Skip to content
Snippets Groups Projects
Commit 8318bd33 authored by Alvin Noe Ladines's avatar Alvin Noe Ladines
Browse files

Fixed conversion

parent c17dbeca
No related branches found
No related tags found
No related merge requests found
......@@ -2,7 +2,6 @@ import os
import numpy as np
import logging
from nomadcore.unit_conversion import unit_conversion
from nomad.datamodel.metainfo.public import section_run, section_system,\
section_single_configuration_calculation, section_method, section_calculation_to_calculation_refs,\
Workflow, Elastic
......@@ -12,10 +11,6 @@ from elasticparser.metainfo.elastic import x_elastic_section_strain_diagrams,\
from elasticparser.elastic_properties import ElasticProperties
a_to_m = unit_conversion.convert_unit_function('angstrom', 'meter')
ha_to_J = unit_conversion.convert_unit(1, 'hartree', 'J')
giga = 10 ** 9
class ElasticParserInterface:
def __init__(self, filepath, archive, logger=None):
......@@ -40,7 +35,7 @@ class ElasticParserInterface:
sec_strain_diagram.x_elastic_strain_diagram_type = 'energy'
sec_strain_diagram.x_elastic_strain_diagram_number_of_eta = len(strain[0])
sec_strain_diagram.x_elastic_strain_diagram_eta_values = strain
sec_strain_diagram.x_elastic_strain_diagram_values = np.array(energy) * ha_to_J
sec_strain_diagram.x_elastic_strain_diagram_values = energy
poly_fit_2 = int((n_strains - 1) / 2)
......@@ -60,8 +55,7 @@ class ElasticParserInterface:
sec_strain_diagram.x_elastic_strain_diagram_polynomial_fit_order = int(fit_order[:-2])
sec_strain_diagram.x_elastic_strain_diagram_number_of_eta = poly_fit.get(fit_order, None)
sec_strain_diagram.x_elastic_strain_diagram_eta_values = energy_fit[diagram_type][0][fit_order]
convert = giga if diagram_type == 'd2e' else ha_to_J
sec_strain_diagram.x_elastic_strain_diagram_values = np.array(energy_fit[diagram_type][1][fit_order]) * convert
sec_strain_diagram.x_elastic_strain_diagram_values = energy_fit[diagram_type][1][fit_order]
elif method == 'stress':
strain, stress = self.properties.get_strain_stress()
......@@ -92,8 +86,7 @@ class ElasticParserInterface:
sec_strain_diagram.x_elastic_strain_diagram_polynomial_fit_order = int(fit_order[:-2])
sec_strain_diagram.x_elastic_strain_diagram_number_of_eta = poly_fit.get(fit_order, None)
sec_strain_diagram.x_elastic_strain_diagram_eta_values = stress_fit[diagram_type][si][0][fit_order]
convert = giga if diagram_type == 'dtn' else ha_to_J
sec_strain_diagram.x_elastic_strain_diagram_values = np.array(stress_fit[diagram_type][si][1][fit_order]) * convert
sec_strain_diagram.x_elastic_strain_diagram_values = np.array(stress_fit[diagram_type][si][1][fit_order])
def parse_elastic_constant(self):
sec_scc = self.archive.section_run[-1].section_single_configuration_calculation[-1]
......@@ -104,26 +97,26 @@ class ElasticParserInterface:
matrices, moduli, eigenvalues = self.properties.get_elastic_constants_order2()
sec_scc.x_elastic_2nd_order_constants_notation_matrix = matrices['voigt']
sec_scc.x_elastic_2nd_order_constants_matrix = matrices['elastic_constant'] * giga
sec_scc.x_elastic_2nd_order_constants_compliance_matrix = matrices['compliance'] / giga
sec_scc.x_elastic_2nd_order_constants_matrix = matrices['elastic_constant']
sec_scc.x_elastic_2nd_order_constants_compliance_matrix = matrices['compliance']
sec_scc.x_elastic_Voigt_bulk_modulus = moduli.get('B_V', moduli.get('K_V')) * giga
sec_scc.x_elastic_Voigt_shear_modulus = moduli['G_V'] * giga
sec_scc.x_elastic_Voigt_bulk_modulus = moduli.get('B_V', moduli.get('K_V'))
sec_scc.x_elastic_Voigt_shear_modulus = moduli['G_V']
sec_scc.x_elastic_Reuss_bulk_modulus = moduli.get('B_R', moduli.get('K_R')) * giga
sec_scc.x_elastic_Reuss_shear_modulus = moduli['G_R'] * giga
sec_scc.x_elastic_Reuss_bulk_modulus = moduli.get('B_R', moduli.get('K_R'))
sec_scc.x_elastic_Reuss_shear_modulus = moduli['G_R']
sec_scc.x_elastic_Hill_bulk_modulus = moduli.get('B_H', moduli.get('K_H')) * giga
sec_scc.x_elastic_Hill_shear_modulus = moduli['G_H'] * giga
sec_scc.x_elastic_Hill_bulk_modulus = moduli.get('B_H', moduli.get('K_H'))
sec_scc.x_elastic_Hill_shear_modulus = moduli['G_H']
sec_scc.x_elastic_Voigt_Young_modulus = moduli['E_V'] * giga
sec_scc.x_elastic_Voigt_Young_modulus = moduli['E_V']
sec_scc.x_elastic_Voigt_Poisson_ratio = moduli['nu_V']
sec_scc.x_elastic_Reuss_Young_modulus = moduli['E_R'] * giga
sec_scc.x_elastic_Reuss_Young_modulus = moduli['E_R']
sec_scc.x_elastic_Reuss_Poisson_ratio = moduli['nu_R']
sec_scc.x_elastic_Hill_Young_modulus = moduli['E_H'] * giga
sec_scc.x_elastic_Hill_Young_modulus = moduli['E_H']
sec_scc.x_elastic_Hill_Poisson_ratio = moduli['nu_H']
sec_scc.x_elastic_eigenvalues = eigenvalues * giga
sec_scc.x_elastic_eigenvalues = eigenvalues
elif order == 3:
elastic_constant = self.properties.get_elastic_constants_order3()
......@@ -141,11 +134,10 @@ class ElasticParserInterface:
symbols, positions, cell = self.properties.get_structure_info()
volume = self.properties.info['equilibrium_volume']
volume = a_to_m(a_to_m(a_to_m(volume)))
sec_system.atom_labels = symbols
sec_system.atom_positions = a_to_m(positions)
sec_system.simulation_cell = a_to_m(cell)
sec_system.atom_positions = positions
sec_system.simulation_cell = cell
sec_system.configuration_periodic_dimensions = [True, True, True]
sec_system.x_elastic_space_group_number = self.properties.info['space_group_number']
sec_system.x_elastic_unit_cell_volume = volume
......
import os
import pint
import numpy as np
from ase import Atoms
......@@ -32,7 +33,7 @@ class ElasticProperties:
Quantity('calculation_method', r'\s*Method of calculation\s*=\s*([-a-zA-Z]+)\s*'),
Quantity('code_name', r'\s*DFT code name\s*=\s*([-a-zA-Z]+)'),
Quantity('space_group_number', r'\s*Space-group number\s*=\s*([0-9]+)'),
Quantity('equilibrium_volume', r'\s*Volume of equilibrium unit cell\s*=\s*([0-9.]+)\s*'),
Quantity('equilibrium_volume', r'\s*Volume of equilibrium unit cell\s*=\s*([0-9.]+)\s*', unit='angstrom ** 3'),
Quantity('max_strain', r'\s*Maximum Lagrangian strain\s*=\s*([0-9.]+)'),
Quantity('n_strains', r'\s*Number of distorted structures\s*=\s*([0-9]+)')
]
......@@ -107,7 +108,12 @@ class ElasticProperties:
structure = Atoms(cell=cellpar, scaled_positions=pos, symbols=sym, pbc=True)
return sym, structure.get_positions(), list(structure.get_cell())
positions = structure.get_positions()
positions = pint.Quantity(positions, 'angstrom')
cell = structure.get_cell()
cell = pint.Quantity(cell, 'angstrom')
return sym, positions, cell
def get_strain_energy(self):
strains, energies = [], []
......@@ -122,6 +128,11 @@ class ElasticProperties:
strains.append(data[0])
energies.append(data[1])
energies = pint.Quantity(energies, 'hartree')
# the peculiarity of the x_elastic_strain_diagram_values metainfo that it does
# not have the energy unit
energies = energies.to('J').magnitude
return strains, energies
def get_strain_stress(self):
......@@ -206,7 +217,7 @@ class ElasticProperties:
eta = {key: val for key, val in eta.items() if val}
val = {key: val for key, val in val.items() if val}
return eta, val
return [eta, val]
def get_energy_fit(self):
energy_fit = dict()
......@@ -216,10 +227,14 @@ class ElasticProperties:
if result is None:
continue
result[1] = {
key: pint.Quantity(val, 'GPa').to('Pa').magnitude for key, val in result[1].items()}
energy_fit['d2e'] = result
result = self._get_fit('Energy-vs-Strain', 'CVe.dat')
if result is not None:
result[1] = {
key: pint.Quantity(val, 'hartree').to('J').magnitude for key, val in result[1].items()}
energy_fit['cross-validation'] = result
return energy_fit
......@@ -232,10 +247,12 @@ class ElasticProperties:
for strain_index in range(1, 7):
result = self._get_fit('Stress-vs-Strain', '%d_dS.dat' % strain_index)
if result is not None:
result[1] = {key: pint.Quantity(val, 'GPa') for key, val in result[1].items()}
stress_fit['dtn'][strain_index - 1] = result
result = self._get_fit('Stress-vs-Strain', '%d_CVe.dat' % strain_index)
if result is not None:
result[1] = {key: pint.Quantity(val, 'hartree') for key, val in result[1].items()}
stress_fit['cross-validation'][strain_index - 1] = result
return stress_fit
......@@ -301,21 +318,25 @@ class ElasticProperties:
str_operation=reshape, dtype=str),
Quantity(
'elastic_constant', r'Elastic constant[\s\S]+in GPa\s*:\s*\n\n([\-\d\.\s\n]+)\n',
str_operation=reshape, dtype=float),
str_operation=reshape, dtype=float, unit='GPa'),
Quantity(
'compliance', r'Elastic compliance[\s\S]+in 1/GPa\s*:\s*\n\n([\-\d\.\s\n]+)\n',
str_operation=reshape, dtype=float)]
str_operation=reshape, dtype=float, unit='1/GPa')]
modulus_names = [
'B_V', 'K_V', 'G_V', 'B_R', 'K_R', 'G_R', 'B_H', 'K_H', 'G_H', 'E_V',
'nu_V', 'E_R', 'nu_R', 'E_H', 'nu_H', 'AVR']
quantities += [
Quantity(name, r'%s\s*=\s*([-+\d+\.]+)\s*' % name) for name in modulus_names]
for name in modulus_names:
unit = None
if name[0:2] in ('B_', 'K_', 'G_', 'E_'):
unit = 'GPa'
quantities.append(
Quantity(name, r'%s\s*=\s*([-+\d+\.]+)\s*' % name, unit=unit))
quantities += [
Quantity(
'eigenvalues', r'Eigenvalues of elastic constant \(stiffness\) matrix:\s*\n+([\-\d\.\n\s]+)\n')]
'eigenvalues', r'Eigenvalues of elastic constant \(stiffness\) matrix:\s*\n+([\-\d\.\n\s]+)\n', unit='GPa')]
parser = UnstructuredTextFileParser(path, quantities)
......
......@@ -330,6 +330,7 @@ class section_system(public.section_system):
x_elastic_unit_cell_volume = Quantity(
type=np.dtype(np.float64),
shape=[],
unit='m ** 3',
description='''
Volume of the equilibrium unit cell
''',
......
0% Loading or .
You are about to add 0 people to the discussion. Proceed with caution.
Please register or to comment