scnc
/
phonopy
mirror of https://github.com/phonopy/phonopy.git
1
0
Fork 0
Code Issues Packages Projects Releases Wiki Activity

Arrange locations of methods in file_IO and interface

This commit is contained in:
Atsushi Togo 2014-11-25 12:18:41 +09:00
parent ca7b168511
commit e8d1d89437
9 changed files with 770 additions and 684 deletions
Show Stats Download Patch File Download Diff File Expand all files Collapse all files

42
anharmonic/file_IO.py
View File

@ -3,8 +3,7 @@ import os
import numpy as np
import h5py
from phonopy.structure.atoms import Atoms
from phonopy.interface import vasp
from phonopy.file_IO import read_force_constant_vasprun_xml, iterparse
from phonopy.interface.vasp import get_forces_from_vasprun_xmls, get_force_constants_from_vasprun_xmls, write_vasp
###########
#
@ -155,7 +154,7 @@ def write_supercells_with_displacements_from_direction_dataset(
positions=positions,
cell=lattice,
pbc=True)
vasp.write_vasp('POSCAR-%05d' % count1, atoms, direct=True)
write_vasp('POSCAR-%05d' % count1, atoms, direct=True)
# YAML
w.write("- number: %5d\n" % (disp1['number'] + 1))
@ -192,7 +191,7 @@ def write_supercells_with_displacements_from_direction_dataset(
cell=lattice,
pbc=True)
if included:
vasp.write_vasp('POSCAR-%05d' % count2, atoms, direct=True)
write_vasp('POSCAR-%05d' % count2, atoms, direct=True)
# YAML
w.write(" - [%20.16f,%20.16f,%20.16f ] # %05d\n" %
@ -261,7 +260,7 @@ def write_supercells_with_three_displacements(supercell,
positions=positions,
cell=lattice,
pbc=True)
vasp.write_vasp('POSCAR-%05d' % count1, atoms, direct=True)
write_vasp('POSCAR-%05d' % count1, atoms, direct=True)
# YAML
w3.write("- number: %5d\n" % (disp1['number'] + 1))
@ -287,7 +286,7 @@ def write_supercells_with_three_displacements(supercell,
positions=positions,
cell=lattice,
pbc=True)
vasp.write_vasp('POSCAR-%05d' % count2, atoms, direct=True)
write_vasp('POSCAR-%05d' % count2, atoms, direct=True)
# YAML
if second_atom_num != disp2['number']:
@ -321,7 +320,7 @@ def write_supercells_with_three_displacements(supercell,
positions=positions,
cell=lattice,
pbc=True)
vasp.write_vasp('POSCAR-%05d' % count3, atoms, direct=True)
write_vasp('POSCAR-%05d' % count3, atoms, direct=True)
# YAML
w4.write(" - [ %20.16f,%20.16f,%20.16f ] # %d \n" %
@ -1268,35 +1267,6 @@ def write_ir_grid_points(mesh,
w.write(" q-point: [ %12.7f, %12.7f, %12.7f ]\n" %
tuple(grid_address[g].astype('double') / mesh))
def get_forces_from_vasprun_xmls(vaspruns, num_atom, index_shift=0):
forces = []
for i, vasprun in enumerate(vaspruns):
print >> sys.stderr, "%d" % (i + 1 + index_shift),
if vasp.is_version528(vasprun):
force_set = vasp.get_forces_vasprun_xml(
iterparse(vasp.VasprunWrapper(vasprun), tag='varray'))
else:
force_set = vasp.get_forces_vasprun_xml(
iterparse(vasprun, tag='varray'))
if force_set.shape[0] == num_atom:
forces.append(force_set)
else:
print "\nNumber of forces in vasprun.xml #%d is wrong." % (i + 1)
sys.exit(1)
print >> sys.stderr
return np.array(forces)
def get_force_constants_from_vasprun_xmls(vasprun_filenames):
force_constants_set = []
for i, filename in enumerate(vasprun_filenames):
print >> sys.stderr, "%d: %s\n" % (i + 1, filename),
force_constants_set.append(
read_force_constant_vasprun_xml(filename)[0])
print >> sys.stderr
return force_constants_set
def parse_yaml(file_yaml):
import yaml
try:

480
phonopy/file_IO.py
View File

@ -36,157 +36,29 @@ import sys
import os
import StringIO
import numpy as np
import phonopy.interface.vasp as vasp
from phonopy.structure.symmetry import Symmetry
from phonopy.harmonic.force_constants import similarity_transformation
from phonopy.structure.atoms import Atoms
from phonopy.cui.settings import fracval
def read_crystal_structure(filename=None,
interface_mode='vasp',
chemical_symbols=None):
if filename is None:
if interface_mode == 'vasp':
unitcell_filename = "POSCAR"
if interface_mode == 'abinit':
unitcell_filename = "unitcell.in"
if interface_mode == 'pwscf':
unitcell_filename = "unitcell.in"
if interface_mode == 'wien2k':
unitcell_filename = "case.struct"
else:
unitcell_filename = filename
if not os.path.exists(unitcell_filename):
if filename is None:
return None, (unitcell_filename + " (default file name)",)
else:
return None, (unitcell_filename,)
#
# FORCE_SETS
#
def write_FORCE_SETS(dataset, filename='FORCE_SETS'):
num_atom = dataset['natom']
displacements = dataset['first_atoms']
forces = [x['forces'] for x in dataset['first_atoms']]
if interface_mode == 'vasp':
from phonopy.interface.vasp import read_vasp
if chemical_symbols is None:
unitcell = read_vasp(unitcell_filename)
else:
unitcell = read_vasp(unitcell_filename, symbols=chemical_symbols)
return unitcell, (unitcell_filename,)
if interface_mode == 'abinit':
from phonopy.interface.abinit import read_abinit
unitcell = read_abinit(unitcell_filename)
return unitcell, (unitcell_filename,)
if interface_mode == 'pwscf':
from phonopy.interface.pwscf import read_pwscf
unitcell, pp_filenames = read_pwscf(unitcell_filename)
return unitcell, (unitcell_filename, pp_filenames)
if interface_mode == 'wien2k':
from phonopy.interface.wien2k import parse_wien2k_struct
unitcell, npts, r0s, rmts = parse_wien2k_struct(unitcell_filename)
return unitcell, (unitcell_filename, npts, r0s, rmts)
def create_FORCE_CONSTANTS(filename, options, log_level):
fc_and_atom_types = read_force_constant_vasprun_xml(filename)
if not fc_and_atom_types:
print
print "\'%s\' dones not contain necessary information." % filename
return 1
# Write FORCE_SETS
fp = open(filename, 'w')
fp.write("%-5d\n" % num_atom)
fp.write("%-5d\n" % len(displacements))
for count, disp in enumerate(displacements):
fp.write("\n%-5d\n" % (disp['number'] + 1))
fp.write("%20.16f %20.16f %20.16f\n" % (tuple(disp['displacement'])))
force_constants, atom_types = fc_and_atom_types
if options.is_hdf5:
try:
import h5py
except ImportError:
print
print "You need to install python-h5py."
return 1
write_force_constants_to_hdf5(force_constants)
if log_level > 0:
print "force_constants.hdf5 has been created from vasprun.xml."
else:
write_FORCE_CONSTANTS(force_constants)
if log_level > 0:
print "FORCE_CONSTANTS has been created from vasprun.xml."
if log_level > 0:
print "Atom types:", atom_types
return 0
def create_FORCE_SETS(interface_mode,
force_filenames,
options,
log_level):
if interface_mode == 'vasp':
displacements = parse_disp_yaml('disp.yaml')
if (interface_mode == 'wien2k' or
interface_mode == 'abinit' or
interface_mode == 'pwscf'):
displacements, supercell = parse_disp_yaml(
'disp.yaml', return_cell=True)
num_disp_files = len(force_filenames)
if options.force_sets_zero_mode:
num_disp_files -= 1
if len(displacements['first_atoms']) != num_disp_files:
print
print ("Number of files to be read don't match "
"to number of displacements in disp.yaml.")
return 1
if interface_mode == 'vasp':
is_created = write_FORCE_SETS_vasp(
force_filenames,
displacements,
filename='FORCE_SETS',
is_zero_point=options.force_sets_zero_mode)
if interface_mode == 'abinit':
print "**********************************************************"
print "**** Abinit FORCE_SETS support is experimental. ****"
print "**** Your feedback would be appreciated. ****"
print "**********************************************************"
is_created = write_FORCE_SETS_abinit(
force_filenames,
displacements,
supercell.get_number_of_atoms(),
filename='FORCE_SETS')
if interface_mode == 'pwscf':
print "**********************************************************"
print "**** Pwscf FORCE_SETS support is experimental. ****"
print "**** Your feedback would be appreciated. ****"
print "**********************************************************"
is_created = write_FORCE_SETS_pwscf(
force_filenames,
displacements,
supercell.get_number_of_atoms(),
filename='FORCE_SETS')
if interface_mode == 'wien2k':
print "**********************************************************"
print "**** Wien2k FORCE_SETS support is experimental. ****"
print "**** Your feedback would be appreciated. ****"
print "**********************************************************"
is_created = write_FORCE_SETS_wien2k(
force_filenames,
displacements,
supercell,
filename='FORCE_SETS',
is_distribute=(not options.is_wien2k_p1),
symprec=options.symprec)
if log_level > 0:
if is_created:
print "FORCE_SETS has been created."
else:
print "FORCE_SETS could not be created."
return 0
def _get_line_ignore_blank(f):
line = f.readline().strip()
if line == '':
line = _get_line_ignore_blank(f)
return line
for f in forces[count]:
fp.write("%15.10f %15.10f %15.10f\n" % (tuple(f)))
def parse_FORCE_SETS(is_translational_invariance=False, filename="FORCE_SETS"):
f = open(filename, 'r')
@ -225,15 +97,13 @@ def _get_set_of_forces(f, is_translational_invariance):
return dataset
def parse_QPOINTS(filename="QPOINTS"):
f = open(filename, 'r')
num_qpoints = int(f.readline().strip())
qpoints = []
for i in range(num_qpoints):
qpoints.append([fracval(x) for x in f.readline().strip().split()])
return np.array(qpoints)
def _get_line_ignore_blank(f):
line = f.readline().strip()
if line == '':
line = _get_line_ignore_blank(f)
return line
def _get_drift_forces(forces, filename=None):
def get_drift_forces(forces, filename=None):
drift_force = np.sum(forces, axis=0) / len(forces)
if filename is None:
print "Drift force"
@ -244,45 +114,7 @@ def _get_drift_forces(forces, filename=None):
return drift_force
def write_FORCE_SETS_abinit(forces_filenames,
displacements,
num_atom,
filename='FORCE_SETS'):
hook = 'cartesian forces (eV/Angstrom)'
for abinit_filename, disp in zip(forces_filenames,
displacements['first_atoms']):
f = open(abinit_filename)
abinit_forces = _collect_forces(f, num_atom, hook, [1, 2, 3])
if not abinit_forces:
return False
drift_force = _get_drift_forces(abinit_forces)
disp['forces'] = np.array(abinit_forces) - drift_force
write_FORCE_SETS(displacements, filename=filename)
return True
def write_FORCE_SETS_pwscf(forces_filenames,
displacements,
num_atom,
filename='FORCE_SETS'):
hook = 'Forces acting on atoms'
for pwscf_filename, disp in zip(forces_filenames,
displacements['first_atoms']):
pwscf_forces = _iter_collect_forces(pwscf_filename,
num_atom,
hook,
[6, 7, 8],
word='force')
drift_force = _get_drift_forces(pwscf_forces)
disp['forces'] = np.array(pwscf_forces) - drift_force
write_FORCE_SETS(displacements, filename=filename)
return True
def _collect_forces(f, num_atom, hook, force_pos, word=None):
def collect_forces(f, num_atom, hook, force_pos, word=None):
for line in f:
if hook in line:
break
@ -308,18 +140,18 @@ def _collect_forces(f, num_atom, hook, force_pos, word=None):
return forces
def _iter_collect_forces(filename,
num_atom,
hook,
force_pos,
word=None,
max_iter=1000):
def iter_collect_forces(filename,
num_atom,
hook,
force_pos,
word=None,
max_iter=1000):
f = open(filename)
forces = []
prev_forces = []
for i in range(max_iter):
forces = _collect_forces(f, num_atom, hook, force_pos, word=word)
forces = collect_forces(f, num_atom, hook, force_pos, word=word)
if not forces:
forces = prev_forces[:]
break
@ -331,140 +163,49 @@ def _iter_collect_forces(filename,
return forces
def write_FORCE_SETS_wien2k(forces_filenames,
displacements,
supercell,
filename='FORCE_SETS',
is_distribute=True,
symprec=1e-5):
import phonopy.interface.wien2k as wien2k
#
# FORCE_CONSTANTS, force_constants.hdf5
#
def write_FORCE_CONSTANTS(force_constants, filename='FORCE_CONSTANTS'):
w = open(filename, 'w')
fc_shape = force_constants.shape
w.write("%4d\n" % (fc_shape[0]))
for i in range(fc_shape[0]):
for j in range(fc_shape[1]):
w.write("%4d%4d\n" % (i+1, j+1))
for vec in force_constants[i][j]:
w.write(("%22.15f"*3 + "\n") % tuple(vec))
w.close()
natom = supercell.get_number_of_atoms()
lattice = supercell.get_cell()
def write_force_constants_to_hdf5(force_constants,
filename='force_constants.hdf5'):
import h5py
w = h5py.File(filename, 'w')
w.create_dataset('force_constants', data=force_constants)
w.close()
for wien2k_filename, disp in zip(forces_filenames,
displacements['first_atoms']):
# Parse wien2k case.scf file
wien2k_forces = wien2k.get_forces_wien2k(wien2k_filename, lattice)
if is_distribute:
force_set = wien2k.distribute_forces(
supercell,
[disp['number'], disp['displacement']],
wien2k_forces,
wien2k_filename,
symprec)
if not force_set:
return False
else:
if not (natom == len(wien2k_forces)):
print "%s contains only forces of %d atoms" % (
wien2k_filename, len(wien2k_forces))
return False
else:
force_set = wien2k_forces
def parse_FORCE_CONSTANTS(filename):
fcfile = open(filename)
num = int((fcfile.readline().strip().split())[0])
force_constants = np.zeros((num, num, 3, 3), dtype=float)
for i in range(num):
for j in range(num):
fcfile.readline()
tensor = []
for k in range(3):
tensor.append([float(x) for x in fcfile.readline().strip().split()])
force_constants[i, j] = np.array(tensor)
drift_force = _get_drift_forces(force_set, filename=wien2k_filename)
disp['forces'] = np.array(force_set) - drift_force
write_FORCE_SETS(displacements, filename=filename)
return True
return force_constants
def read_force_constants_hdf5(filename="force_constants.hdf5"):
import h5py
f = h5py.File(filename)
return f[f.keys()[0]][:]
def iterparse(fname, tag=None):
try:
import xml.etree.cElementTree as etree
for event, elem in etree.iterparse(fname):
if tag is None or elem.tag == tag:
yield event, elem
except ImportError:
print "Python 2.5 or later is needed."
print "For creating FORCE_SETS file with Python 2.4, you can use",
print "phonopy 1.8.5.1 with python-lxml ."
sys.exit(1)
def write_FORCE_SETS_vasp(forces_filenames,
displacements,
filename='FORCE_SETS',
is_zero_point=False,
verbose=True):
if verbose:
print "counter (file index):",
num_atom = displacements['natom']
count = 0
are_files_correct = True
if is_zero_point:
force_files = forces_filenames[1:]
if vasp.is_version528(forces_filenames[0]):
zero_forces = vasp.get_forces_vasprun_xml(iterparse(
vasp.VasprunWrapper(forces_filenames[0]), tag='varray'))
else:
zero_forces = vasp.get_forces_vasprun_xml(
iterparse(forces_filenames[0], tag='varray'))
if verbose:
print "%d" % (count + 1),
count += 1
if not _check_forces(zero_forces, num_atom, forces_filenames[0]):
are_files_correct = False
else:
force_files = forces_filenames
zero_forces = None
for i, disp in enumerate(displacements['first_atoms']):
if vasp.is_version528(force_files[i]):
disp['forces'] = vasp.get_forces_vasprun_xml(iterparse(
vasp.VasprunWrapper(force_files[i]), tag='varray'))
else:
disp['forces'] = vasp.get_forces_vasprun_xml(
iterparse(force_files[i], tag='varray'))
if verbose:
print "%d" % (count + 1),
count += 1
if not _check_forces(disp['forces'], num_atom, force_files[i]):
are_files_correct = False
if verbose:
print
write_FORCE_SETS(displacements,
filename=filename,
zero_forces=zero_forces)
return are_files_correct
def _check_forces(forces, num_atom, filename):
if len(forces) != num_atom:
print " \"%s\" does not contain necessary information." % filename,
return False
else:
return True
def write_FORCE_SETS(dataset, filename='FORCE_SETS', zero_forces=None):
num_atom = dataset['natom']
displacements = dataset['first_atoms']
forces = [x['forces'] for x in dataset['first_atoms']]
# Write FORCE_SETS
fp = open(filename, 'w')
fp.write("%-5d\n" % num_atom)
fp.write("%-5d\n" % len(displacements))
for count, disp in enumerate(displacements):
fp.write("\n%-5d\n" % (disp['number'] + 1))
fp.write("%20.16f %20.16f %20.16f\n" % (tuple(disp['displacement'])))
# Subtract residual forces
if zero_forces is not None:
forces[count] -= zero_forces
for f in forces[count]:
fp.write("%15.10f %15.10f %15.10f\n" % (tuple(f)))
#
# disp.yaml
#
def parse_disp_yaml(filename="disp.yaml", return_cell=False):
try:
import yaml
@ -507,17 +248,6 @@ def parse_disp_yaml(filename="disp.yaml", return_cell=False):
else:
return new_dataset
def parse_DISP(filename='DISP'):
disp = open(filename)
displacements = []
for line in disp:
if line.strip() != '':
a = line.split()
displacements.append(
[int(a[0])-1, float(a[1]), float(a[2]), float(a[3])])
return displacements
# Write disp.yaml
def write_disp_yaml(displacements, supercell, directions=None,
filename='disp.yaml'):
file = open(filename, 'w')
@ -543,57 +273,33 @@ def write_disp_yaml(displacements, supercell, directions=None,
(v[0], v[1], v[2]))
file.close()
# Write FORCE_CONSTANTS
def write_FORCE_CONSTANTS(force_constants, filename='FORCE_CONSTANTS'):
w = open(filename, 'w')
fc_shape = force_constants.shape
w.write("%4d\n" % (fc_shape[0]))
for i in range(fc_shape[0]):
for j in range(fc_shape[1]):
w.write("%4d%4d\n" % (i+1, j+1))
for vec in force_constants[i][j]:
w.write(("%22.15f"*3 + "\n") % tuple(vec))
w.close()
#
# DISP (old phonopy displacement format)
#
def parse_DISP(filename='DISP'):
disp = open(filename)
displacements = []
for line in disp:
if line.strip() != '':
a = line.split()
displacements.append(
[int(a[0])-1, float(a[1]), float(a[2]), float(a[3])])
return displacements
def write_force_constants_to_hdf5(force_constants,
filename='force_constants.hdf5'):
import h5py
w = h5py.File(filename, 'w')
w.create_dataset('force_constants', data=force_constants)
w.close()
#
# QPOINTS
#
def parse_QPOINTS(filename="QPOINTS"):
f = open(filename, 'r')
num_qpoints = int(f.readline().strip())
qpoints = []
for i in range(num_qpoints):
qpoints.append([fracval(x) for x in f.readline().strip().split()])
return np.array(qpoints)
# Read FORCE_CONSTANTS
def parse_FORCE_CONSTANTS(filename):
fcfile = open(filename)
num = int((fcfile.readline().strip().split())[0])
force_constants = np.zeros((num, num, 3, 3), dtype=float)
for i in range(num):
for j in range(num):
fcfile.readline()
tensor = []
for k in range(3):
tensor.append([float(x) for x in fcfile.readline().strip().split()])
force_constants[i, j] = np.array(tensor)
return force_constants
def read_force_constant_vasprun_xml(filename):
if vasp.is_version528(filename):
vasprun = iterparse(vasp.VasprunWrapper(filename))
else:
vasprun = iterparse(filename)
return vasp.get_force_constants_vasprun_xml(vasprun)
def read_force_constant_OUTCAR(filename):
return vasp.get_force_constants_OUTCAR(filename)
def read_force_constants_hdf5(filename="force_constants.hdf5"):
import h5py
f = h5py.File(filename)
return f[f.keys()[0]][:]
# Read BORN
#
# BORN
#
def parse_BORN(primitive, symprec=1e-5, is_symmetry=True, filename="BORN"):
f = open(filename, 'r')
symmetry = Symmetry(primitive, symprec=symprec, is_symmetry=is_symmetry)

156
phonopy/interface/__init__.py
View File

@ -0,0 +1,156 @@
# Copyright (C) 2014 Atsushi Togo
# All rights reserved.
#
# This file is part of phonopy.
#
# Redistribution and use in source and binary forms, with or without
# modification, are permitted provided that the following conditions
# are met:
#
# * Redistributions of source code must retain the above copyright
# notice, this list of conditions and the following disclaimer.
#
# * Redistributions in binary form must reproduce the above copyright
# notice, this list of conditions and the following disclaimer in
# the documentation and/or other materials provided with the
# distribution.
#
# * Neither the name of the phonopy project nor the names of its
# contributors may be used to endorse or promote products derived
# from this software without specific prior written permission.
#
# THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS
# "AS IS" AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT
# LIMITED TO, THE IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS
# FOR A PARTICULAR PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL THE
# COPYRIGHT HOLDER OR CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT,
# INCIDENTAL, SPECIAL, EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING,
# BUT NOT LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES;
# LOSS OF USE, DATA, OR PROFITS; OR BUSINESS INTERRUPTION) HOWEVER
# CAUSED AND ON ANY THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT
# LIABILITY, OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN
# ANY WAY OUT OF THE USE OF THIS SOFTWARE, EVEN IF ADVISED OF THE
# POSSIBILITY OF SUCH DAMAGE.
import os
from phonopy.file_IO import parse_disp_yaml, write_FORCE_SETS
def read_crystal_structure(filename=None,
interface_mode='vasp',
chemical_symbols=None):
if filename is None:
if interface_mode == 'vasp':
unitcell_filename = "POSCAR"
if interface_mode == 'abinit':
unitcell_filename = "unitcell.in"
if interface_mode == 'pwscf':
unitcell_filename = "unitcell.in"
if interface_mode == 'wien2k':
unitcell_filename = "case.struct"
else:
unitcell_filename = filename
if not os.path.exists(unitcell_filename):
if filename is None:
return None, (unitcell_filename + " (default file name)",)
else:
return None, (unitcell_filename,)
if interface_mode == 'vasp':
from phonopy.interface.vasp import read_vasp
if chemical_symbols is None:
unitcell = read_vasp(unitcell_filename)
else:
unitcell = read_vasp(unitcell_filename, symbols=chemical_symbols)
return unitcell, (unitcell_filename,)
if interface_mode == 'abinit':
from phonopy.interface.abinit import read_abinit
unitcell = read_abinit(unitcell_filename)
return unitcell, (unitcell_filename,)
if interface_mode == 'pwscf':
from phonopy.interface.pwscf import read_pwscf
unitcell, pp_filenames = read_pwscf(unitcell_filename)
return unitcell, (unitcell_filename, pp_filenames)
if interface_mode == 'wien2k':
from phonopy.interface.wien2k import parse_wien2k_struct
unitcell, npts, r0s, rmts = parse_wien2k_struct(unitcell_filename)
return unitcell, (unitcell_filename, npts, r0s, rmts)
def create_FORCE_SETS(interface_mode,
force_filenames,
options,
log_level):
if interface_mode == 'vasp':
displacements = parse_disp_yaml(filename='disp.yaml')
if (interface_mode == 'wien2k' or
interface_mode == 'abinit' or
interface_mode == 'pwscf'):
displacements, supercell = parse_disp_yaml(filename='disp.yaml',
return_cell=True)
num_disp_files = len(force_filenames)
if options.force_sets_zero_mode:
num_disp_files -= 1
if len(displacements['first_atoms']) != num_disp_files:
print
print ("Number of files to be read don't match "
"to number of displacements in disp.yaml.")
return 1
if interface_mode == 'vasp':
from phonopy.interface.vasp import parse_set_of_forces
is_parsed = parse_set_of_forces(
displacements,
force_filenames,
is_zero_point=options.force_sets_zero_mode)
if interface_mode == 'abinit':
from phonopy.interface.abinit import parse_set_of_forces
print "**********************************************************"
print "**** Abinit FORCE_SETS support is experimental. ****"
print "**** Your feedback would be appreciated. ****"
print "**********************************************************"
is_parsed = parse_set_of_forces(
displacements,
force_filenames,
supercell.get_number_of_atoms())
if interface_mode == 'pwscf':
from phonopy.interface.pwscf import parse_set_of_forces
print "**********************************************************"
print "**** Pwscf FORCE_SETS support is experimental. ****"
print "**** Your feedback would be appreciated. ****"
print "**********************************************************"
is_parsed = parse_set_of_forces(
displacements,
force_filenames,
supercell.get_number_of_atoms())
if interface_mode == 'wien2k':
from phonopy.interface.wien2k import parse_set_of_forces
print "**********************************************************"
print "**** Wien2k FORCE_SETS support is experimental. ****"
print "**** Your feedback would be appreciated. ****"
print "**********************************************************"
is_parsed = parse_set_of_forces(
displacements,
force_filenames,
supercell,
is_distribute=(not options.is_wien2k_p1),
symprec=options.symprec)
if is_parsed:
write_FORCE_SETS(displacements,
filename='FORCE_SETS')
if log_level > 0:
if is_parsed:
print "FORCE_SETS has been created."
else:
print "FORCE_SETS could not be created."
return 0

17
phonopy/interface/abinit.py
View File

@ -35,11 +35,28 @@
import sys
import numpy as np
from phonopy.file_IO import collect_forces, get_drift_forces
from phonopy.interface.vasp import get_scaled_positions_lines
from phonopy.units import Bohr
from phonopy.cui.settings import fracval
from phonopy.structure.atoms import Atoms
def parse_set_of_forces(displacements,
forces_filenames,
num_atom):
hook = 'cartesian forces (eV/Angstrom)'
for abinit_filename, disp in zip(forces_filenames,
displacements['first_atoms']):
f = open(abinit_filename)
abinit_forces = collect_forces(f, num_atom, hook, [1, 2, 3])
if not abinit_forces:
return False
drift_force = get_drift_forces(abinit_forces)
disp['forces'] = np.array(abinit_forces) - drift_force
return True
def read_abinit(filename):
abinit_in = AbinitIn(open(filename).readlines())
tags = abinit_in.get_variables()

17
phonopy/interface/pwscf.py
View File

@ -35,11 +35,28 @@
import sys
import numpy as np
from phonopy.file_IO import iter_collect_forces, get_drift_forces
from phonopy.interface.vasp import get_scaled_positions_lines
from phonopy.units import Bohr
from phonopy.cui.settings import fracval
from phonopy.structure.atoms import Atoms, symbol_map
def parse_set_of_forces(displacements,
forces_filenames,
num_atom):
hook = 'Forces acting on atoms'
for pwscf_filename, disp in zip(forces_filenames,
displacements['first_atoms']):
pwscf_forces = iter_collect_forces(pwscf_filename,
num_atom,
hook,
[6, 7, 8],
word='force')
drift_force = get_drift_forces(pwscf_forces)
disp['forces'] = np.array(pwscf_forces) - drift_force
return True
def read_pwscf(filename):
pwscf_in = PwscfIn(open(filename).readlines())
tags = pwscf_in.get_tags()

631
phonopy/interface/vasp.py
View File

@ -39,29 +39,190 @@ from phonopy.structure.atoms import Atoms, symbol_map, atom_data
from phonopy.structure.cells import get_primitive, get_supercell
from phonopy.structure.symmetry import Symmetry
from phonopy.harmonic.force_constants import similarity_transformation
from phonopy.file_IO import write_FORCE_SETS, write_force_constants_to_hdf5, write_FORCE_CONSTANTS
class VasprunWrapper(object):
"""VasprunWrapper class
This is used to avoid VASP 5.2.8 vasprun.xml defect at PRECFOCK,
xml parser stops with error.
"""
def __init__(self, filename):
self.f = open(filename)
def parse_set_of_forces(displacements,
forces_filenames,
is_zero_point=False,
verbose=True):
def read(self, size=None):
element = self.f.next()
if element.find("PRECFOCK") == -1:
return element
if verbose:
print "counter (file index):",
num_atom = displacements['natom']
count = 0
is_parsed = True
if is_zero_point:
force_files = forces_filenames[1:]
if _is_version528(forces_filenames[0]):
zero_forces = _get_forces_vasprun_xml(_iterparse(
VasprunWrapper(forces_filenames[0]), tag='varray'))
else:
return "<i type=\"string\" name=\"PRECFOCK\"></i>"
zero_forces = _get_forces_vasprun_xml(
_iterparse(forces_filenames[0], tag='varray'))
def is_version528(filename):
for line in open(filename):
if '\"version\"' in line:
if '5.2.8' in line:
return True
else:
return False
if verbose:
print "%d" % (count + 1),
count += 1
if not _check_forces(zero_forces, num_atom, forces_filenames[0]):
is_parsed = False
else:
force_files = forces_filenames
for i, disp in enumerate(displacements['first_atoms']):
if _is_version528(force_files[i]):
disp['forces'] = _get_forces_vasprun_xml(_iterparse(
VasprunWrapper(force_files[i]), tag='varray'))
else:
disp['forces'] = _get_forces_vasprun_xml(
_iterparse(force_files[i], tag='varray'))
if is_zero_point:
disp['forces'] -= zero_forces
if verbose:
print "%d" % (count + 1),
count += 1
if not _check_forces(disp['forces'], num_atom, force_files[i]):
is_parsed = False
if verbose:
print
return is_parsed
def _check_forces(forces, num_atom, filename):
if len(forces) != num_atom:
print " \"%s\" does not contain necessary information." % filename,
return False
else:
return True
def create_FORCE_CONSTANTS(filename, options, log_level):
fc_and_atom_types = read_force_constant_vasprun_xml(filename)
if not fc_and_atom_types:
print
print "\'%s\' dones not contain necessary information." % filename
return 1
force_constants, atom_types = fc_and_atom_types
if options.is_hdf5:
try:
import h5py
except ImportError:
print
print "You need to install python-h5py."
return 1
write_force_constants_to_hdf5(force_constants)
if log_level > 0:
print "force_constants.hdf5 has been created from vasprun.xml."
else:
write_FORCE_CONSTANTS(force_constants)
if log_level > 0:
print "FORCE_CONSTANTS has been created from vasprun.xml."
if log_level > 0:
print "Atom types:", atom_types
return 0
#
# read VASP POSCAR
#
def read_vasp(filename, symbols=None):
lines = open(filename).readlines()
return _get_atoms_from_poscar(lines, symbols)
def read_vasp_from_strings(strings, symbols=None):
return _get_atoms_from_poscar(
StringIO.StringIO(strings).readlines(), symbols)
def _get_atoms_from_poscar(lines, symbols):
line1 = [x for x in lines[0].split()]
if _is_exist_symbols(line1):
symbols = line1
scale = float(lines[1])
cell = []
for i in range(2, 5):
cell.append([float(x) for x in lines[i].split()[:3]])
cell = np.array(cell) * scale
try:
num_atoms = np.array([int(x) for x in lines[5].split()])
line_at = 6
except ValueError:
symbols = [x for x in lines[5].split()]
num_atoms = np.array([int(x) for x in lines[6].split()])
line_at = 7
expaned_symbols = _expand_symbols(num_atoms, symbols)
if lines[line_at][0].lower() == 's':
line_at += 1
is_scaled = True
if (lines[line_at][0].lower() == 'c' or
lines[line_at][0].lower() == 'k'):
is_scaled = False
line_at += 1
positions = []
for i in range(line_at, line_at + num_atoms.sum()):
positions.append([float(x) for x in lines[i].split()[:3]])
if is_scaled:
atoms = Atoms(symbols=expaned_symbols,
cell=cell,
scaled_positions=positions)
else:
atoms = Atoms(symbols=expaned_symbols,
cell=cell,
positions=positions)
return atoms
def _is_exist_symbols(symbols):
for s in symbols:
if not (s in symbol_map):
return False
return True
def _expand_symbols(num_atoms, symbols=None):
expanded_symbols = []
is_symbols = True
if symbols is None:
is_symbols = False
else:
if len(symbols) != len(num_atoms):
is_symbols = False
else:
for s in symbols:
if not s in symbol_map:
is_symbols = False
break
if is_symbols:
for s, num in zip(symbols, num_atoms):
expanded_symbols += [s] * num
else:
for i, num in enumerate(num_atoms):
expanded_symbols += [atom_data[i+1][1]] * num
return expanded_symbols
#
# write vasp POSCAR
#
def write_vasp(filename, atoms, direct=True):
lines = _get_vasp_structure(atoms, direct=direct)
f = open(filename, 'w')
f.write(lines)
def write_supercells_with_displacements(supercell,
cells_with_displacements):
@ -69,10 +230,9 @@ def write_supercells_with_displacements(supercell,
for i, cell in enumerate(cells_with_displacements):
write_vasp('POSCAR-%03d' % (i + 1), cell, direct=True)
write_magnetic_moments(supercell)
_write_magnetic_moments(supercell)
def write_magnetic_moments(cell):
def _write_magnetic_moments(cell):
magmoms = cell.get_magnetic_moments()
if magmoms is not None:
w = open("MAGMOM", 'w')
@ -85,97 +245,62 @@ def write_magnetic_moments(cell):
w.write("\n")
w.close()
def get_scaled_positions_lines(scaled_positions):
lines = ""
for i, vec in enumerate(scaled_positions):
for x in (vec - np.rint(vec)):
if float('%20.16f' % x) < 0.0:
lines += "%20.16f" % (x + 1.0)
else:
lines += "%20.16f" % (x)
if i < len(scaled_positions) - 1:
lines += "\n"
def get_forces_vasprun_xml(vasprun):
"""
vasprun = etree.iterparse(filename, tag='varray')
"""
forces = []
for event, element in vasprun:
if element.attrib['name'] == 'forces':
for v in element:
forces.append([float(x) for x in v.text.split()])
return np.array(forces)
return lines
def get_force_constants_vasprun_xml(vasprun):
fc_tmp = None
num_atom = 0
for event, element in vasprun:
if num_atom == 0:
(atom_types,
masses,
num_atom) = get_atom_types_from_vasprun_xml(element)
def _get_vasp_structure(atoms, direct=True):
(num_atoms,
symbols,
scaled_positions,
sort_list) = _sort_positions_by_symbols(atoms.get_chemical_symbols(),
atoms.get_scaled_positions())
lines = ""
for s in symbols:
lines += "%s " % s
lines += "\n"
lines += " 1.0\n"
for a in atoms.get_cell():
lines += " %22.16f%22.16f%22.16f\n" % tuple(a)
lines += ("%4d" * len(num_atoms)) % tuple(num_atoms)
lines += "\n"
lines += "Direct\n"
lines += get_scaled_positions_lines(scaled_positions)
# Get Hessian matrix (normalized by masses)
if element.tag == 'varray':
if element.attrib['name'] == 'hessian':
fc_tmp = []
for v in element.findall('./v'):
fc_tmp.append([float(x) for x in v.text.strip().split()])
if fc_tmp is None:
return False
else:
fc_tmp = np.array(fc_tmp)
if fc_tmp.shape != (num_atom * 3, num_atom * 3):
return False
# num_atom = fc_tmp.shape[0] / 3
force_constants = np.zeros((num_atom, num_atom, 3, 3), dtype='double')
return lines
for i in range(num_atom):
for j in range(num_atom):
force_constants[i, j] = fc_tmp[i*3:(i+1)*3, j*3:(j+1)*3]
# Inverse normalization by atomic weights
for i in range(num_atom):
for j in range(num_atom):
force_constants[i, j] *= -np.sqrt(masses[i] * masses[j])
def _get_reduced_symbols(symbols):
reduced_symbols = []
for s in symbols:
if not (s in reduced_symbols):
reduced_symbols.append(s)
return reduced_symbols
return force_constants, atom_types
def _sort_positions_by_symbols(symbols, positions):
reduced_symbols = _get_reduced_symbols(symbols)
sorted_positions = []
sort_list = []
num_atoms = np.zeros(len(reduced_symbols), dtype=int)
for i, rs in enumerate(reduced_symbols):
for j, (s, p) in enumerate(zip(symbols, positions)):
if rs == s:
sorted_positions.append(p)
sort_list.append(j)
num_atoms[i] += 1
return num_atoms, reduced_symbols, np.array(sorted_positions), sort_list
def get_atom_types_from_vasprun_xml(element):
atom_types = []
masses = []
num_atom = 0
if element.tag == 'array':
if 'name' in element.attrib:
if element.attrib['name'] == 'atomtypes':
for rc in element.findall('./set/rc'):
atom_info = [x.text for x in rc.findall('./c')]
num_atom += int(atom_info[0])
atom_types.append(atom_info[1].strip())
masses += ([float(atom_info[2])] * int(atom_info[0]))
return atom_types, masses, num_atom
def get_force_constants_OUTCAR(filename):
file = open(filename)
while 1:
line = file.readline()
if line == '':
print "Force constants could not be found."
return 0
if line[:19] == " SECOND DERIVATIVES":
break
file.readline()
num_atom = int(((file.readline().split())[-1].strip())[:-1])
fc_tmp = []
for i in range(num_atom * 3):
fc_tmp.append([float(x) for x in (file.readline().split())[1:]])
fc_tmp = np.array(fc_tmp)
force_constants = np.zeros((num_atom, num_atom, 3, 3), dtype=float)
for i in range(num_atom):
for j in range(num_atom):
force_constants[i, j] = -fc_tmp[i*3:(i+1)*3, j*3:(j+1)*3]
return force_constants
#
# Non-analytical term
#
def get_born_OUTCAR(poscar_filename="POSCAR",
outcar_filename="OUTCAR",
primitive_axis=np.eye(3),
@ -291,155 +416,177 @@ def _symmetrize_tensor(tensor, symmetry_operations):
return sum_tensor / len(symmetry_operations)
#
# read VASP POSCAR
# vasprun.xml handling
#
def read_vasp(filename, symbols=None):
lines = open(filename).readlines()
return _get_atoms_from_poscar(lines, symbols)
class VasprunWrapper(object):
"""VasprunWrapper class
This is used to avoid VASP 5.2.8 vasprun.xml defect at PRECFOCK,
xml parser stops with error.
"""
def __init__(self, filename):
self.f = open(filename)
def read_vasp_from_strings(strings, symbols=None):
return _get_atoms_from_poscar(
StringIO.StringIO(strings).readlines(), symbols)
def _get_atoms_from_poscar(lines, symbols):
line1 = [x for x in lines[0].split()]
if _is_exist_symbols(line1):
symbols = line1
scale = float(lines[1])
cell = []
for i in range(2, 5):
cell.append([float(x) for x in lines[i].split()[:3]])
cell = np.array(cell) * scale
try:
num_atoms = np.array([int(x) for x in lines[5].split()])
line_at = 6
except ValueError:
symbols = [x for x in lines[5].split()]
num_atoms = np.array([int(x) for x in lines[6].split()])
line_at = 7
expaned_symbols = _expand_symbols(num_atoms, symbols)
if lines[line_at][0].lower() == 's':
line_at += 1
is_scaled = True
if (lines[line_at][0].lower() == 'c' or
lines[line_at][0].lower() == 'k'):
is_scaled = False
line_at += 1
positions = []
for i in range(line_at, line_at + num_atoms.sum()):
positions.append([float(x) for x in lines[i].split()[:3]])
if is_scaled:
atoms = Atoms(symbols=expaned_symbols,
cell=cell,
scaled_positions=positions)
else:
atoms = Atoms(symbols=expaned_symbols,
cell=cell,
positions=positions)
return atoms
def _is_exist_symbols(symbols):
for s in symbols:
if not (s in symbol_map):
return False
return True
def _expand_symbols(num_atoms, symbols=None):
expanded_symbols = []
is_symbols = True
if symbols is None:
is_symbols = False
else:
if len(symbols) != len(num_atoms):
is_symbols = False
def read(self, size=None):
element = self.f.next()
if element.find("PRECFOCK") == -1:
return element
else:
for s in symbols:
if not s in symbol_map:
is_symbols = False
break
if is_symbols:
for s, num in zip(symbols, num_atoms):
expanded_symbols += [s] * num
else:
for i, num in enumerate(num_atoms):
expanded_symbols += [atom_data[i+1][1]] * num
return "<i type=\"string\" name=\"PRECFOCK\"></i>"
return expanded_symbols
def _iterparse(fname, tag=None):
try:
import xml.etree.cElementTree as etree
for event, elem in etree.iterparse(fname):
if tag is None or elem.tag == tag:
yield event, elem
except ImportError:
print "Python 2.5 or later is needed."
print "For creating FORCE_SETS file with Python 2.4, you can use",
print "phonopy 1.8.5.1 with python-lxml ."
sys.exit(1)
#
# write vasp POSCAR
#
def get_scaled_positions_lines(scaled_positions):
lines = ""
for i, vec in enumerate(scaled_positions):
for x in (vec - np.rint(vec)):
if float('%20.16f' % x) < 0.0:
lines += "%20.16f" % (x + 1.0)
def _is_version528(filename):
for line in open(filename):
if '\"version\"' in line:
if '5.2.8' in line:
return True
else:
lines += "%20.16f" % (x)
if i < len(scaled_positions) - 1:
lines += "\n"
return False
return lines
def read_force_constant_vasprun_xml(filename):
if _is_version528(filename):
vasprun = _iterparse(VasprunWrapper(filename))
else:
vasprun = _iterparse(filename)
return get_force_constants_vasprun_xml(vasprun)
def get_vasp_structure(atoms, direct=True):
(num_atoms,
symbols,
scaled_positions,
sort_list) = _sort_positions_by_symbols(atoms.get_chemical_symbols(),
atoms.get_scaled_positions())
lines = ""
for s in symbols:
lines += "%s " % s
lines += "\n"
lines += " 1.0\n"
for a in atoms.get_cell():
lines += " %22.16f%22.16f%22.16f\n" % tuple(a)
lines += ("%4d" * len(num_atoms)) % tuple(num_atoms)
lines += "\n"
lines += "Direct\n"
lines += get_scaled_positions_lines(scaled_positions)
def get_force_constants_vasprun_xml(vasprun):
fc_tmp = None
num_atom = 0
for event, element in vasprun:
if num_atom == 0:
(atom_types,
masses,
num_atom) = _get_atom_types_from_vasprun_xml(element)
return lines
# Get Hessian matrix (normalized by masses)
if element.tag == 'varray':
if element.attrib['name'] == 'hessian':
fc_tmp = []
for v in element.findall('./v'):
fc_tmp.append([float(x) for x in v.text.strip().split()])
if fc_tmp is None:
return False
else:
fc_tmp = np.array(fc_tmp)
if fc_tmp.shape != (num_atom * 3, num_atom * 3):
return False
# num_atom = fc_tmp.shape[0] / 3
force_constants = np.zeros((num_atom, num_atom, 3, 3), dtype='double')
def write_vasp(filename, atoms, direct=True):
lines = get_vasp_structure(atoms, direct=direct)
f = open(filename, 'w')
f.write(lines)
for i in range(num_atom):
for j in range(num_atom):
force_constants[i, j] = fc_tmp[i*3:(i+1)*3, j*3:(j+1)*3]
# Inverse normalization by atomic weights
for i in range(num_atom):
for j in range(num_atom):
force_constants[i, j] *= -np.sqrt(masses[i] * masses[j])
def _get_reduced_symbols(symbols):
reduced_symbols = []
for s in symbols:
if not (s in reduced_symbols):
reduced_symbols.append(s)
return reduced_symbols
return force_constants, atom_types
def get_forces_from_vasprun_xmls(vaspruns, num_atom, index_shift=0):
forces = []
for i, vasprun in enumerate(vaspruns):
print >> sys.stderr, "%d" % (i + 1 + index_shift),
def _sort_positions_by_symbols(symbols, positions):
reduced_symbols = _get_reduced_symbols(symbols)
sorted_positions = []
sort_list = []
num_atoms = np.zeros(len(reduced_symbols), dtype=int)
for i, rs in enumerate(reduced_symbols):
for j, (s, p) in enumerate(zip(symbols, positions)):
if rs == s:
sorted_positions.append(p)
sort_list.append(j)
num_atoms[i] += 1
return num_atoms, reduced_symbols, np.array(sorted_positions), sort_list
if _is_version528(vasprun):
force_set = _get_forces_vasprun_xml(
_iterparse(VasprunWrapper(vasprun), tag='varray'))
else:
force_set = _get_forces_vasprun_xml(
_iterparse(vasprun, tag='varray'))
if force_set.shape[0] == num_atom:
forces.append(force_set)
else:
print "\nNumber of forces in vasprun.xml #%d is wrong." % (i + 1)
sys.exit(1)
print >> sys.stderr
return np.array(forces)
def get_force_constants_from_vasprun_xmls(vasprun_filenames):
force_constants_set = []
for i, filename in enumerate(vasprun_filenames):
print >> sys.stderr, "%d: %s\n" % (i + 1, filename),
force_constants_set.append(
read_force_constant_vasprun_xml(filename)[0])
print >> sys.stderr
return force_constants_set
def _get_forces_vasprun_xml(vasprun):
"""
vasprun = etree.iterparse(filename, tag='varray')
"""
forces = []
for event, element in vasprun:
if element.attrib['name'] == 'forces':
for v in element:
forces.append([float(x) for x in v.text.split()])
return np.array(forces)
def _get_atom_types_from_vasprun_xml(element):
atom_types = []
masses = []
num_atom = 0
if element.tag == 'array':
if 'name' in element.attrib:
if element.attrib['name'] == 'atomtypes':
for rc in element.findall('./set/rc'):
atom_info = [x.text for x in rc.findall('./c')]
num_atom += int(atom_info[0])
atom_types.append(atom_info[1].strip())
masses += ([float(atom_info[2])] * int(atom_info[0]))
return atom_types, masses, num_atom
#
# OUTCAR handling (obsolete)
#
def read_force_constant_OUTCAR(filename):
return get_force_constants_OUTCAR(filename)
def get_force_constants_OUTCAR(filename):
file = open(filename)
while 1:
line = file.readline()
if line == '':
print "Force constants could not be found."
return 0
if line[:19] == " SECOND DERIVATIVES":
break
file.readline()
num_atom = int(((file.readline().split())[-1].strip())[:-1])
fc_tmp = []
for i in range(num_atom * 3):
fc_tmp.append([float(x) for x in (file.readline().split())[1:]])
fc_tmp = np.array(fc_tmp)
force_constants = np.zeros((num_atom, num_atom, 3, 3), dtype=float)
for i in range(num_atom):
for j in range(num_atom):
force_constants[i, j] = -fc_tmp[i*3:(i+1)*3, j*3:(j+1)*3]
return force_constants
if __name__ == '__main__':
import sys
atoms = read_vasp(sys.argv[1])
write_vasp('%s-new' % sys.argv[1], atoms)

71
phonopy/interface/wien2k.py
View File

@ -33,11 +33,82 @@
# POSSIBILITY OF SUCH DAMAGE.
import numpy as np
from phonopy.file_IO import write_FORCE_SETS, get_drift_forces
from phonopy.structure.atoms import Atoms
from phonopy.structure.symmetry import Symmetry
from phonopy.structure.cells import get_angles, get_cell_parameters
from phonopy.harmonic.force_constants import similarity_transformation
def parse_set_of_forces(displacements,
forces_filenames,
supercell,
is_distribute=True,
symprec=1e-5):
natom = supercell.get_number_of_atoms()
lattice = supercell.get_cell()
for wien2k_filename, disp in zip(forces_filenames,
displacements['first_atoms']):
# Parse wien2k case.scf file
wien2k_forces = get_forces_wien2k(wien2k_filename, lattice)
if is_distribute:
force_set = distribute_forces(
supercell,
[disp['number'], disp['displacement']],
wien2k_forces,
wien2k_filename,
symprec)
if not force_set:
return False
else:
if not (natom == len(wien2k_forces)):
print "%s contains only forces of %d atoms" % (
wien2k_filename, len(wien2k_forces))
return False
else:
force_set = wien2k_forces
drift_force = get_drift_forces(force_set, filename=wien2k_filename)
disp['forces'] = np.array(force_set) - drift_force
return True
def create_FORCE_SETS(forces_filenames,
displacements,
supercell,
filename='FORCE_SETS',
is_distribute=True,
symprec=1e-5):
natom = supercell.get_number_of_atoms()
lattice = supercell.get_cell()
for wien2k_filename, disp in zip(forces_filenames,
displacements['first_atoms']):
# Parse wien2k case.scf file
wien2k_forces = get_forces_wien2k(wien2k_filename, lattice)
if is_distribute:
force_set = distribute_forces(
supercell,
[disp['number'], disp['displacement']],
wien2k_forces,
wien2k_filename,
symprec)
if not force_set:
return False
else:
if not (natom == len(wien2k_forces)):
print "%s contains only forces of %d atoms" % (
wien2k_filename, len(wien2k_forces))
return False
else:
force_set = wien2k_forces
drift_force = get_drift_forces(force_set, filename=wien2k_filename)
disp['forces'] = np.array(force_set) - drift_force
write_FORCE_SETS(displacements, filename=filename)
return True
def get_wien2k_struct(cell, npts, r0s, rmts):
num_atom = cell.get_number_of_atoms()

2
scripts/dispmanager
View File

@ -42,7 +42,7 @@
import sys, os
import numpy as np
from phonopy.structure.atoms import Atoms
import phonopy.interface.vasp as vasp
from phonopy.interface.vasp import write_vasp
import phonopy.file_IO as file_IO
try:

38
scripts/phonopy
View File

@ -45,6 +45,7 @@ from phonopy.cui.show_symmetry import check_symmetry
from phonopy.units import *
from phonopy.version import phonopy_version
from phonopy.structure.cells import print_cell, determinant
from phonopy.interface import create_FORCE_SETS, read_crystal_structure
# AA is created at http://www.network-science.de/ascii/.
@ -465,19 +466,20 @@ if options.loglevel is not None:
#
if options.wien2k_mode:
interface_mode = 'wien2k'
import phonopy.interface.wien2k as wien2k
from phonopy.interface.wien2k import write_supercells_with_displacements
factor = Wien2kToTHz
elif options.abinit_mode:
interface_mode = 'abinit'
import phonopy.interface.abinit as abinit
from phonopy.interface.abinit import write_supercells_with_displacements
factor = AbinitToTHz
elif options.pwscf_mode:
interface_mode = 'pwscf'
import phonopy.interface.pwscf as pwscf
from phonopy.interface.pwscf import write_supercells_with_displacements
factor = PwscfToTHz
else:
interface_mode = 'vasp'
import phonopy.interface.vasp as vasp
from phonopy.interface.vasp import write_supercells_with_displacements
from phonopy.interface.vasp import create_FORCE_CONSTANTS
factor = VaspToTHz
# Show title
@ -490,10 +492,10 @@ if options.force_sets_mode or options.force_sets_zero_mode:
file_exists('disp.yaml', log_level)
for filename in args:
file_exists(filename, log_level)
error_num = file_IO.create_FORCE_SETS(interface_mode,
args,
options,
log_level)
error_num = create_FORCE_SETS(interface_mode,
args,
options,
log_level)
if log_level > 0:
print_end()
sys.exit(error_num)
@ -502,9 +504,9 @@ if options.force_sets_mode or options.force_sets_zero_mode:
if options.force_constants_mode:
if len(args) > 0:
file_exists(args[0], log_level)
error_num = file_IO.create_FORCE_CONSTANTS(args[0],
options,
log_level)
error_num = create_FORCE_CONSTANTS(args[0],
options,
log_level)
else:
print_error_message("Please specify vasprun.xml.")
error_num = 1
@ -533,7 +535,7 @@ if options.is_graph_save:
# Initialization #
##################
unitcell, optional_structure_file_information = file_IO.read_crystal_structure(
unitcell, optional_structure_file_information = read_crystal_structure(
filename=settings.get_cell_filename(),
interface_mode=interface_mode,
chemical_symbols=settings.get_chemical_symbols())
@ -713,7 +715,7 @@ if run_mode == 'displacements':
if interface_mode == 'wien2k':
npts, r0s, rmts = optional_structure_file_information[1:4]
wien2k.write_supercells_with_displacements(
write_supercells_with_displacements(
supercell,
cells_with_disps,
npts,
@ -722,14 +724,14 @@ if run_mode == 'displacements':
settings.get_supercell_matrix(),
filename=unitcell_filename)
elif interface_mode == 'abinit':
abinit.write_supercells_with_displacements(supercell, cells_with_disps)
write_supercells_with_displacements(supercell, cells_with_disps)
elif interface_mode == 'pwscf':
pp_filenames = optional_structure_file_information[1]
pwscf.write_supercells_with_displacements(supercell,
cells_with_disps,
pp_filenames)
write_supercells_with_displacements(supercell,
cells_with_disps,
pp_filenames)
else:
vasp.write_supercells_with_displacements(supercell, cells_with_disps)
write_supercells_with_displacements(supercell, cells_with_disps)
if log_level > 0:
print