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

Minor update for QHA

This commit is contained in:
Atsushi Togo 2018-08-07 14:57:37 +09:00
parent ee4075a2fa
commit e2cb25c480
1 changed files with 23 additions and 32 deletions
Show Stats Download Patch File Download Diff File Expand all files Collapse all files

55
scripts/phonopy-qha
View File

@ -49,7 +49,7 @@ def get_options():
is_graph_plot=False,
is_graph_save=False,
is_bulk_modulus_only=False,
read_efe_file=False,
efe_file=None,
eos="vinet",
thin_number=10,
tmax=1000.0)
@ -69,9 +69,9 @@ def get_options():
"--pressure", dest="pressure", type=float,
help="Pressure in GPa")
parser.add_argument(
"--read-efe", "--read-electronic-free-energy",
dest="read_efe_file", action="store_true",
help="Read electronic free energies from file instead of v-e data")
"--efe", "--electronic-free-energy",
dest="efe_file", nargs=1,
help="Read electronic free energies at temperatures and volumes")
parser.add_argument(
"-s", "--save", dest="is_graph_save", action="store_true",
help="Save plot data in pdf")
@ -105,8 +105,8 @@ def main(args):
# Show bulk modulus of v-e data
if args.is_bulk_modulus_only:
if args.read_efe_file:
print("--read-efe optin can't be used with -b option")
if args.efe_file:
print("--efe optin can't be used with -b option")
sys.exit(1)
volumes, electronic_energies = read_v_e(args.filenames[0])
@ -126,34 +126,30 @@ def main(args):
########################
# Read data from files #
########################
if args.read_efe_file:
_temperatures, electronic_energies = read_efe(args.filenames[0])
else:
volumes, electronic_energies = read_v_e(args.filenames[0])
electronic_energies += volumes * args.pressure / EVAngstromToGPa
# Check number of files in e-v.dat case
if not args.read_efe_file and len(volumes) != len(args.filenames[1:]):
print("The number of thermal_properites.yaml files (%d) "
"is inconsisten with" % len(args.filenames[1:]))
print("the number of e-v data (%d)." % len(volumes))
volumes, electronic_energies = read_v_e(args.filenames[0])
electronic_energies += volumes * args.pressure / EVAngstromToGPa
# Check number of files in e-v.dat case
if len(volumes) != len(args.filenames[1:]):
print("The number of thermal_properites.yaml files (%d) "
"is inconsisten with" % len(args.filenames[1:]))
print("the number of e-v data (%d)." % len(volumes))
sys.exit(1)
if args.efe_file:
_temperatures, electronic_energies = read_efe(args.efe_file[0])
if len(volumes) != electronic_energies.shape[1]:
print("%s and %s are inconsistent for the volume points."
% (args.filenames[0], args.efe_file[0]))
sys.exit(1)
(temperatures,
cv,
entropy,
fe_phonon,
_volumes,
num_modes,
num_integrated_modes) = read_thermal_properties_yaml(args.filenames[1:])
if args.read_efe_file:
if len(_volumes) == electronic_energies.shape[1]:
volumes = np.array(_volumes)
else:
print("Unit cell volumes could not be collected from "
"thermal_properties.yaml files. Use the latset phonopy-qha "
"and re-calculate thermal_properties.yaml's")
sys.exit(1)
if args.efe_file:
if ((len(temperatures) >= len(_temperatures) and
(np.abs(temperatures[:len(_temperatures)] - _temperatures)
> 1e-5).any()) or
@ -162,12 +158,7 @@ def main(args):
> 1e-5).any())):
print("Inconsistency is found in temperatures in %s and "
"thermal_properties.yaml files." % args.filenames[0])
else:
if len(_volumes) == len(volumes):
if (np.abs(volumes - _volumes) > 1e-4).any():
print("Inconsistency is found in unti cell volumes in %s "
"and thermal_properties.yaml files" % args.filenames[0])
sys.exit(1)
sys.exit(1)
########################################################
# Treatment of thermal properties with imaginary modes #
@ -181,7 +172,7 @@ def main(args):
else:
indices = range(len(volumes))
if args.read_efe_file:
if args.efe_file:
electronic_energies = electronic_energies[:, indices]
else:
electronic_energies = electronic_energies[indices]