From 49681648d8e198e5f3d887088a6ef20d1a89737a Mon Sep 17 00:00:00 2001 From: "Jose M. Pizarro" <112697669+JosePizarro3@users.noreply.github.com> Date: Tue, 28 May 2024 14:53:46 +0200 Subject: [PATCH] Add DielectricFunction properties (#72) * Fix variables module and added Frequency Commented out Variables.points.shape to fix when points are refs to a specific NumericalSettings section (like KMesh(Variables) or KLinePath(Variables)) Deleted refs to numerical settings under KMesh(Variables) and KLinePath(Variables) * Changed to plural in the list of properties under Outputs Added AbsorptionSpectrum to the list * Added Permittivity property Added extract_absorption_spectra from the imaginary part of the permittivity Added AbsorptionSpectrum.axis to differentiate the principal axis to which is defined from * Added get_variables utils function Changed code to use this utils function * Added testing for Permittivity and get_variables Fixed testing --- src/nomad_simulations/numerical_settings.py | 2 +- src/nomad_simulations/outputs.py | 20 ++- src/nomad_simulations/properties/__init__.py | 4 +- .../properties/permittivity.py | 118 +++++++++++++ .../properties/spectral_profile.py | 93 ++++++---- src/nomad_simulations/utils/__init__.py | 7 +- src/nomad_simulations/utils/utils.py | 25 ++- src/nomad_simulations/variables.py | 52 +++--- tests/conftest.py | 6 +- tests/test_numerical_settings.py | 2 - tests/test_outputs.py | 16 +- tests/test_permittivity.py | 167 ++++++++++++++++++ tests/test_spectral_profile.py | 74 ++++---- tests/test_utils.py | 34 +++- 14 files changed, 493 insertions(+), 127 deletions(-) create mode 100644 src/nomad_simulations/properties/permittivity.py create mode 100644 tests/test_permittivity.py diff --git a/src/nomad_simulations/numerical_settings.py b/src/nomad_simulations/numerical_settings.py index 163c0897..e2a2f5eb 100644 --- a/src/nomad_simulations/numerical_settings.py +++ b/src/nomad_simulations/numerical_settings.py @@ -336,7 +336,7 @@ class KMesh(Mesh): description=""" Dictionary containing the high-symmetry point labels and their values in units of `reciprocal_lattice_vectors`. E.g., in a cubic lattice: - high_symmetry_points ={ + high_symmetry_points = { 'Gamma': [0, 0, 0], 'X': [0.5, 0, 0], 'Y': [0, 0.5, 0], diff --git a/src/nomad_simulations/outputs.py b/src/nomad_simulations/outputs.py index ac1ca674..46436a6f 100644 --- a/src/nomad_simulations/outputs.py +++ b/src/nomad_simulations/outputs.py @@ -33,7 +33,9 @@ HoppingMatrix, ElectronicBandGap, ElectronicDensityOfStates, - XASSpectra, + AbsorptionSpectrum, + XASSpectrum, + Permittivity, ) @@ -63,23 +65,27 @@ class Outputs(ArchiveSection): # List of properties # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # - fermi_level = SubSection(sub_section=FermiLevel.m_def, repeats=True) + fermi_levels = SubSection(sub_section=FermiLevel.m_def, repeats=True) - chemical_potential = SubSection(sub_section=ChemicalPotential.m_def, repeats=True) + chemical_potentials = SubSection(sub_section=ChemicalPotential.m_def, repeats=True) - crystal_field_splitting = SubSection( + crystal_field_splittings = SubSection( sub_section=CrystalFieldSplitting.m_def, repeats=True ) - hopping_matrix = SubSection(sub_section=HoppingMatrix.m_def, repeats=True) + hopping_matrices = SubSection(sub_section=HoppingMatrix.m_def, repeats=True) - electronic_band_gap = SubSection(sub_section=ElectronicBandGap.m_def, repeats=True) + electronic_band_gaps = SubSection(sub_section=ElectronicBandGap.m_def, repeats=True) electronic_dos = SubSection( sub_section=ElectronicDensityOfStates.m_def, repeats=True ) - xas_spectra = SubSection(sub_section=XASSpectra.m_def, repeats=True) + permittivities = SubSection(sub_section=Permittivity.m_def, repeats=True) + + absorption_spectra = SubSection(sub_section=AbsorptionSpectrum.m_def, repeats=True) + + xas_spectra = SubSection(sub_section=XASSpectrum.m_def, repeats=True) # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # diff --git a/src/nomad_simulations/properties/__init__.py b/src/nomad_simulations/properties/__init__.py index e529ed38..889196cd 100644 --- a/src/nomad_simulations/properties/__init__.py +++ b/src/nomad_simulations/properties/__init__.py @@ -22,6 +22,8 @@ SpectralProfile, DOSProfile, ElectronicDensityOfStates, - XASSpectra, + AbsorptionSpectrum, + XASSpectrum, ) from .hopping_matrix import HoppingMatrix, CrystalFieldSplitting +from .permittivity import Permittivity diff --git a/src/nomad_simulations/properties/permittivity.py b/src/nomad_simulations/properties/permittivity.py new file mode 100644 index 00000000..f4e4511c --- /dev/null +++ b/src/nomad_simulations/properties/permittivity.py @@ -0,0 +1,118 @@ +# +# Copyright The NOMAD Authors. +# +# This file is part of NOMAD. See https://nomad-lab.eu for further info. +# +# Licensed under the Apache License, Version 2.0 (the "License"); +# you may not use this file except in compliance with the License. +# You may obtain a copy of the License at +# +# http://www.apache.org/licenses/LICENSE-2.0 +# +# Unless required by applicable law or agreed to in writing, software +# distributed under the License is distributed on an "AS IS" BASIS, +# WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied. +# See the License for the specific language governing permissions and +# limitations under the License. +# + +import numpy as np +from structlog.stdlib import BoundLogger +from typing import Optional, List + +from nomad.metainfo import Quantity, Section, Context, MEnum + +from nomad_simulations.physical_property import PhysicalProperty +from nomad_simulations.utils import get_variables +from nomad_simulations.variables import Frequency, KMesh +from nomad_simulations.properties.spectral_profile import AbsorptionSpectrum + + +# TODO add `DielectricStrength` when we have examples and understand how to extract it from the `Permittivity` tensor. + + +class Permittivity(PhysicalProperty): + """ + Response of the material to polarize in the presence of an electric field. + + Alternative names: `DielectricFunction`. + """ + + iri = 'http://fairmat-nfdi.eu/taxonomy/Permittivity' + + type = Quantity( + type=MEnum('static', 'dynamic'), + description=""" + Type of permittivity which allows to identify if the permittivity depends on the frequency or not. + """, + ) + + value = Quantity( + type=np.complex128, + # unit='joule', # TODO check units (they have to match `SpectralProfile.value`) + description=""" + Value of the permittivity tensor. If the value does not depend on the scattering vector `q`, then we + can extract the optical absorption spectrum from the imaginary part of the permittivity tensor (this is also called + macroscopic dielectric function). + """, + ) + + # ? We need use cases to understand if we need to define contributions to the permittivity tensor. + # ? `ionic` and `electronic` contributions are common in the literature. + + def __init__( + self, m_def: Section = None, m_context: Context = None, **kwargs + ) -> None: + super().__init__(m_def, m_context, **kwargs) + self.rank = [3, 3] + self.name = self.m_def.name + self._axes_map = ['xx', 'yy', 'zz'] + + def resolve_type(self) -> str: + frequencies = get_variables(self.variables, Frequency) + if len(frequencies) > 0: + return 'dynamic' + return 'static' + + def extract_absorption_spectra( + self, logger: BoundLogger + ) -> Optional[List[AbsorptionSpectrum]]: + """ + Extract the absorption spectrum from the imaginary part of the permittivity tensor. + """ + # If the `pemittivity` depends on the scattering vector `q`, then we cannot extract the absorption spectrum + q_mesh = get_variables(self.variables, KMesh) + if len(q_mesh) > 0: + logger.warning( + 'The `permittivity` depends on the scattering vector `q`, so that we cannot extract the absorption spectrum.' + ) + return None + # Extract the `Frequency` variable to extract the absorption spectrum + frequencies = get_variables(self.variables, Frequency) + if len(frequencies) == 0: + logger.warning( + 'The `permittivity` does not have a `Frequency` variable to extract the absorption spectrum.' + ) + return None + # Define the `absorption_spectra` for each principal direction along the diagonal of the `Permittivity.value` as the imaginary part + spectra = [] + for i in range(3): + val = self.value[:, i, i].imag + absorption_spectrum = AbsorptionSpectrum( + axis=self._axes_map[i], variables=frequencies + ) + absorption_spectrum.value = val + absorption_spectrum.physical_property_ref = self + spectra.append(absorption_spectrum) + return spectra + + def normalize(self, archive, logger) -> None: + super().normalize(archive, logger) + + # Resolve the `type` of permittivity + self.type = self.resolve_type() + + # `AbsorptionSpectrum` extraction + absorption_spectra = self.extract_absorption_spectra(logger) + if absorption_spectra is not None: + self.m_parent.absorption_spectrum = absorption_spectra diff --git a/src/nomad_simulations/properties/spectral_profile.py b/src/nomad_simulations/properties/spectral_profile.py index 4b57ba60..03bd60ac 100644 --- a/src/nomad_simulations/properties/spectral_profile.py +++ b/src/nomad_simulations/properties/spectral_profile.py @@ -22,9 +22,9 @@ import pint from nomad import config -from nomad.metainfo import Quantity, SubSection, Section, Context +from nomad.metainfo import Quantity, SubSection, Section, Context, MEnum -from ..utils import get_sibling_section +from ..utils import get_sibling_section, get_variables from ..physical_property import PhysicalProperty from ..variables import Energy2 as Energy from ..atoms_state import AtomsState, OrbitalsState @@ -49,24 +49,6 @@ def __init__( super().__init__(m_def, m_context, **kwargs) self.rank = [] - def _get_energy_points(self, logger: BoundLogger) -> Optional[pint.Quantity]: - """ - Gets the `points` of the `Energy` variable if the required `Energy` variable is present in the `variables`. - - Args: - logger (BoundLogger): The logger to log messages. - - Returns: - (Optional[pint.Quantity]): The `points` of the `Energy` variable. - """ - for var in self.variables: - if isinstance(var, Energy): - return var.points - logger.error( - 'The required `Energy` variable is not present in the `variables`.' - ) - return None - def is_valid_spectral_profile(self) -> bool: """ Check if the spectral profile is valid, i.e., if all `value` are defined positive. @@ -450,6 +432,14 @@ def generate_from_projected_dos( if self.projected_dos is None or len(self.projected_dos) == 0: return None + # Extract `Energy` variables + energies = get_variables(self.variables, Energy) + if len(energies) != 1: + logger.warning( + 'The `ElectronicDensityOfStates` does not contain an `Energy` variable to extract the DOS.' + ) + return None + # We distinguish between orbital and atom `projected_dos` orbital_projected = self.extract_projected_dos('orbital', logger) atom_projected = self.extract_projected_dos('atom', logger) @@ -468,7 +458,7 @@ def generate_from_projected_dos( atom_dos = DOSProfile( name=f'atom {ref.chemical_symbol}', entity_ref=ref, - variables=data[0].variables, + variables=energies, ) orbital_values = [ dos.value.magnitude for dos in data @@ -493,9 +483,10 @@ def normalize(self, archive, logger) -> None: super().normalize(archive, logger) # Initial check to see if `variables` contains the required `Energy` variable - energies = self._get_energy_points(logger) - if energies is None: + energies = get_variables(self.variables, Energy) + if len(energies) != 0: return + energies = energies[0].points # Resolve `fermi_level` from a sibling section with respect to `ElectronicDensityOfStates` fermi_level = get_sibling_section( @@ -528,25 +519,49 @@ def normalize(self, archive, logger) -> None: self.value = value_from_pdos -class XASSpectra(SpectralProfile): +class AbsorptionSpectrum(SpectralProfile): + """ """ + + # ! implement `iri` and `rank` as part of `m_def = Section()` + + axis = Quantity( + type=MEnum('xx', 'yy', 'zz'), + description=""" + Axis of the absorption spectrum. This is related with the polarization direction, and can be seen as the + principal term in the tensor `Permittivity.value` (see permittivity.py module). + """, + ) + + def __init__( + self, m_def: Section = None, m_context: Context = None, **kwargs + ) -> None: + super().__init__(m_def, m_context, **kwargs) + # Set the name of the section + self.name = self.m_def.name + + def normalize(self, archive, logger) -> None: + super().normalize(archive, logger) + + +class XASSpectrum(AbsorptionSpectrum): """ - X-ray Absorption Spectra (XAS). + X-ray Absorption Spectrum (XAS). """ # ! implement `iri` and `rank` as part of `m_def = Section()` - xanes_spectra = SubSection( - sub_section=SpectralProfile.m_def, + xanes_spectrum = SubSection( + sub_section=AbsorptionSpectrum.m_def, description=""" - X-ray Absorption Near Edge Structure (XANES) spectra. + X-ray Absorption Near Edge Structure (XANES) spectrum. """, repeats=False, ) - exafs_spectra = SubSection( - sub_section=SpectralProfile.m_def, + exafs_spectrum = SubSection( + sub_section=AbsorptionSpectrum.m_def, description=""" - Extended X-ray Absorption Fine Structure (EXAFS) spectra. + Extended X-ray Absorption Fine Structure (EXAFS) spectrum. """, repeats=False, ) @@ -560,22 +575,24 @@ def __init__( def generate_from_contributions(self, logger: BoundLogger) -> None: """ - Generate the `value` of the XAS spectra by concatenating the XANES and EXAFS contributions. It also concatenates + Generate the `value` of the XAS spectrum by concatenating the XANES and EXAFS contributions. It also concatenates the `Energy` grid points of the XANES and EXAFS parts. Args: logger (BoundLogger): The logger to log messages. """ # TODO check if this method is general enough - if self.xanes_spectra is not None and self.exafs_spectra is not None: + if self.xanes_spectrum is not None and self.exafs_spectrum is not None: # Concatenate XANE and EXAFS `Energy` grid points - xanes_energies = self.xanes_spectra._get_energy_points(logger) - exafs_energies = self.exafs_spectra._get_energy_points(logger) - if xanes_energies is None or exafs_energies is None: + xanes_variables = get_variables(self.xanes_spectrum.variables, Energy) + exafs_variables = get_variables(self.exafs_spectrum.variables, Energy) + if len(xanes_variables) == 0 or len(exafs_variables) == 0: logger.warning( 'Could not extract the `Energy` grid points from XANES or EXAFS.' ) return + xanes_energies = xanes_variables[0].points + exafs_energies = exafs_variables[0].points if xanes_energies.max() > exafs_energies.min(): logger.warning( 'The XANES `Energy` grid points are not below the EXAFS `Energy` grid points.' @@ -587,7 +604,7 @@ def generate_from_contributions(self, logger: BoundLogger) -> None: # Concatenate XANES and EXAFS `value` if they have the same shape ['n_energies'] try: self.value = np.concatenate( - [self.xanes_spectra.value, self.exafs_spectra.value] + [self.xanes_spectrum.value, self.exafs_spectrum.value] ) except ValueError: logger.warning( @@ -599,6 +616,6 @@ def normalize(self, archive, logger) -> None: if self.value is None: logger.info( - 'The `XASSpectra.value` is not stored. We will attempt to obtain it by combining the XANES and EXAFS parts if these are present.' + 'The `XASSpectrum.value` is not stored. We will attempt to obtain it by combining the XANES and EXAFS parts if these are present.' ) self.generate_from_contributions(logger) diff --git a/src/nomad_simulations/utils/__init__.py b/src/nomad_simulations/utils/__init__.py index a0a83520..0c96821f 100644 --- a/src/nomad_simulations/utils/__init__.py +++ b/src/nomad_simulations/utils/__init__.py @@ -16,4 +16,9 @@ # See the License for the specific language governing permissions and # limitations under the License. -from .utils import get_sibling_section, RussellSaundersState, is_not_representative +from .utils import ( + get_sibling_section, + RussellSaundersState, + is_not_representative, + get_variables, +) diff --git a/src/nomad_simulations/utils/utils.py b/src/nomad_simulations/utils/utils.py index e2fe582f..1048d921 100644 --- a/src/nomad_simulations/utils/utils.py +++ b/src/nomad_simulations/utils/utils.py @@ -18,7 +18,7 @@ # from math import factorial -from typing import Optional +from typing import Optional, List from structlog.stdlib import BoundLogger from nomad.datamodel.data import ArchiveSection @@ -129,3 +129,26 @@ def is_not_representative(model_system, logger: BoundLogger): logger.warning('The `ModelSystem` was not found to be representative.') return True return False + + +# cannot define typing with `Variables` due to circular import issue +def get_variables( + variables: Optional[List[ArchiveSection]], variable_cls: ArchiveSection +) -> List[ArchiveSection]: + """ + Get the list of variables which are of type `variable_cls` and appear under `variables`. + + Args: + variables (Optional[List[Variables]]): The list of variables to check. + variable_cls (Variables): The class of the variables to get. + + Returns: + (List[Variables]): The list of variables which are of type `variable_cls`. + """ + if variables is None: + return [] + result = [] + for var in variables: + if isinstance(var, variable_cls): + result.append(var) + return result diff --git a/src/nomad_simulations/variables.py b/src/nomad_simulations/variables.py index 0fe97d75..d3fd1859 100644 --- a/src/nomad_simulations/variables.py +++ b/src/nomad_simulations/variables.py @@ -52,7 +52,7 @@ class Variables(ArchiveSection): points = Quantity( type=np.float64, - shape=['n_points'], + # shape=['n_points'], # ! if defined, this breaks using `points` as refs (e.g., `KMesh.points`) description=""" Points in which the variable is discretized. It might be overwritten with specific units. """, @@ -157,20 +157,34 @@ def normalize(self, archive, logger) -> None: super().normalize(archive, logger) +class Frequency(Variables): + """ """ + + points = Quantity( + type=np.float64, + unit='joule', + shape=['n_points'], + description=""" + Points in which the frequency is discretized in joules. + """, + ) + + def __init__( + self, m_def: Section = None, m_context: Context = None, **kwargs + ) -> None: + super().__init__(m_def, m_context, **kwargs) + self.name = self.m_def.name + + def normalize(self, archive, logger) -> None: + super().normalize(archive, logger) + + class KMesh(Variables): """ K-point mesh over which the physical property is calculated. This is used to define `ElectronicEigenvalues(PhysicalProperty)` and other k-space properties. The `points` are obtained from a refernece to the `NumericalSettings` section, `KMesh(NumericalSettings)`. """ - k_mesh_settings_ref = Quantity( - type=KMeshSettings, - description=""" - Reference to the `KMesh(NumericalSettings)` section in the `ModelMethod` section. This reference is useful - to extract `points` and, then, obtain the shape of `value` of the `PhysicalProperty`. - """, - ) - points = Quantity( type=KMeshSettings.points, description=""" @@ -185,26 +199,12 @@ def __init__( self.name = self.m_def.name def normalize(self, archive, logger) -> None: - # Extracting `points` from the `k_mesh_settings_ref` BEFORE doing `super().normalize()` - if self.k_mesh_settings_ref is None: - logger.error('`k_mesh_settings_ref` is not defined.') - return - self.points = self.k_mesh_settings_ref # ref to `points` - super().normalize(archive, logger) class KLinePath(Variables): """ """ - k_line_path_settings_ref = Quantity( - type=KLinePathSettings, - description=""" - Reference to the `KLinePath(NumericalSettings)` section in the `ModelMethod.KMesh` section. This reference is useful - to extract `points` and, then, obtain the shape of `value` of the `PhysicalProperty`. - """, - ) - points = Quantity( type=KLinePathSettings.points, description=""" @@ -219,10 +219,4 @@ def __init__( self.name = self.m_def.name def normalize(self, archive, logger) -> None: - # Extracting `points` from the `k_line_path_settings_ref` BEFORE doing `super().normalize()` - if self.k_line_path_settings_ref is None: - logger.error('`k_line_path_settings_ref` is not defined.') - return - self.points = self.k_line_path_settings_ref # ref to `points` - super().normalize(archive, logger) diff --git a/tests/conftest.py b/tests/conftest.py index 39107094..2f62bb31 100644 --- a/tests/conftest.py +++ b/tests/conftest.py @@ -170,16 +170,16 @@ def generate_scf_electronic_band_gap_template( for i in range(1, n_scf_steps): value = 1 + sum([1 / (10**j) for j in range(1, i + 1)]) scf_step = Outputs( - electronic_band_gap=[ElectronicBandGap(value=value * ureg.joule)] + electronic_band_gaps=[ElectronicBandGap(value=value * ureg.joule)] ) scf_outputs.scf_steps.append(scf_step) # Add a SCF calculated PhysicalProperty - scf_outputs.electronic_band_gap.append(ElectronicBandGap(value=value * ureg.joule)) + scf_outputs.electronic_band_gaps.append(ElectronicBandGap(value=value * ureg.joule)) # and a `SelfConsistency` ref section scf_params = SelfConsistency( threshold_change=threshold_change, threshold_change_unit='joule' ) - scf_outputs.electronic_band_gap[0].self_consistency_ref = scf_params + scf_outputs.electronic_band_gaps[0].self_consistency_ref = scf_params return scf_outputs diff --git a/tests/test_numerical_settings.py b/tests/test_numerical_settings.py index 3b780da7..11d94ba4 100644 --- a/tests/test_numerical_settings.py +++ b/tests/test_numerical_settings.py @@ -252,7 +252,6 @@ def test_resolve_k_line_density( reciprocal_lattice_vectors = k_space.reciprocal_lattice_vectors k_mesh = k_space.k_mesh[0] model_systems = simulation.model_system - # Applying method `get_k_line_density` get_k_line_density_value = k_mesh.get_k_line_density( reciprocal_lattice_vectors=reciprocal_lattice_vectors, logger=logger @@ -264,7 +263,6 @@ def test_resolve_k_line_density( ) else: assert get_k_line_density_value == result_get_k_line_density - # Applying method `resolve_k_line_density` k_line_density = k_mesh.resolve_k_line_density( model_systems=model_systems, diff --git a/tests/test_outputs.py b/tests/test_outputs.py index aac05ced..d4efaaa2 100644 --- a/tests/test_outputs.py +++ b/tests/test_outputs.py @@ -44,9 +44,9 @@ def test_is_scf_converged(self, threshold_change: float, result: bool): threshold_change=threshold_change ) is_scf_converged = scf_outputs.resolve_is_scf_converged( - property_name='electronic_band_gap', + property_name='electronic_band_gaps', i_property=0, - phys_property=scf_outputs.electronic_band_gap[0], + phys_property=scf_outputs.electronic_band_gaps[0], logger=logger, ) assert is_scf_converged == result @@ -60,18 +60,18 @@ def test_extract_spin_polarized_properties(self): # No spin-polarized band gap band_gap_non_spin_polarized = ElectronicBandGap(variables=[]) band_gap_non_spin_polarized.value = 2.0 * ureg.joule - outputs.electronic_band_gap.append(band_gap_non_spin_polarized) - band_gaps = outputs.extract_spin_polarized_property('electronic_band_gap') + outputs.electronic_band_gaps.append(band_gap_non_spin_polarized) + band_gaps = outputs.extract_spin_polarized_property('electronic_band_gaps') assert band_gaps == [] # Spin-polarized band gaps band_gap_spin_1 = ElectronicBandGap(variables=[], spin_channel=0) band_gap_spin_1.value = 1.0 * ureg.joule - outputs.electronic_band_gap.append(band_gap_spin_1) + outputs.electronic_band_gaps.append(band_gap_spin_1) band_gap_spin_2 = ElectronicBandGap(variables=[], spin_channel=1) band_gap_spin_2.value = 1.5 * ureg.joule - outputs.electronic_band_gap.append(band_gap_spin_2) - band_gaps = outputs.extract_spin_polarized_property('electronic_band_gap') + outputs.electronic_band_gaps.append(band_gap_spin_2) + band_gaps = outputs.extract_spin_polarized_property('electronic_band_gaps') assert len(band_gaps) == 2 assert band_gaps[0].value.magnitude == 1.0 assert band_gaps[1].value.magnitude == 1.5 @@ -89,4 +89,4 @@ def test_normalize(self, threshold_change: float, result: bool): ) scf_outputs.normalize(EntryArchive(), logger) - assert scf_outputs.electronic_band_gap[0].is_scf_converged == result + assert scf_outputs.electronic_band_gaps[0].is_scf_converged == result diff --git a/tests/test_permittivity.py b/tests/test_permittivity.py new file mode 100644 index 00000000..8e85d10c --- /dev/null +++ b/tests/test_permittivity.py @@ -0,0 +1,167 @@ +# +# Copyright The NOMAD Authors. +# +# This file is part of NOMAD. See https://nomad-lab.eu for further info. +# +# Licensed under the Apache License, Version 2.0 (the "License"); +# you may not use this file except in compliance with the License. +# You may obtain a copy of the License at +# +# http://www.apache.org/licenses/LICENSE-2.0 +# +# Unless required by applicable law or agreed to in writing, software +# distributed under the License is distributed on an "AS IS" BASIS, +# WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied. +# See the License for the specific language governing permissions and +# limitations under the License. +# + +import pytest +import numpy as np +from typing import List, Optional + +from nomad.units import ureg +from nomad.datamodel import EntryArchive + +from . import logger + +from nomad_simulations import Simulation +from nomad_simulations.model_system import ModelSystem, AtomicCell +from nomad_simulations.atoms_state import AtomsState +from nomad_simulations.outputs import Outputs +from nomad_simulations.properties import Permittivity +from nomad_simulations.variables import Variables, KMesh, Frequency + +from .conftest import generate_k_space_simulation + + +class TestPermittivity: + """ + Test the `Permittivity` class defined in `properties/permittivity.py`. + """ + + # ! Include this initial `test_default_quantities` method when testing your PhysicalProperty classes + def test_default_quantities(self): + """ + Test the default quantities assigned when creating an instance of the `Permittivity` class. + """ + permittivity = Permittivity() + assert permittivity.iri == 'http://fairmat-nfdi.eu/taxonomy/Permittivity' + assert permittivity.name == 'Permittivity' + assert permittivity.rank == [3, 3] + + @pytest.mark.parametrize( + 'kmesh_grid, variables, result', + [ + (None, None, 'static'), + ([], None, 'static'), + ([3, 3, 3], None, 'static'), + ([3, 3, 3], Frequency(), 'dynamic'), + ], + ) + def test_resolve_type( + self, + kmesh_grid: Optional[List[int]], + variables: Optional[Variables], + result: str, + ): + """ + Test the `resolve_type` method. + """ + permittivity = Permittivity() + + # Generating `KMesh` numerical settings section + simulation = generate_k_space_simulation(grid=kmesh_grid) + k_mesh_settings = simulation.model_method[0].numerical_settings[0].k_mesh[0] + k_mesh_settings.label = 'q-mesh' + k_mesh_settings.center = 'Monkhorst-Pack' + k_mesh_settings.points, _ = k_mesh_settings.resolve_points_and_offset(logger) + if kmesh_grid is not None and len(kmesh_grid) > 0: + kmesh_variables = KMesh(points=k_mesh_settings) + permittivity.variables.append(kmesh_variables) + + # Adding variables to `permittivity` property + if variables is not None and len(variables) > 0: + permittivity.variables.append(variables) + assert permittivity.resolve_type() == result + + @pytest.mark.parametrize( + 'kmesh_grid, frequencies_points, value, result', + [ + # Empty case + (None, None, None, None), + # No `variables` + ([], [], np.eye(3) * (1 + 1j), None), + # If `KMesh` is defined we cannot extract absorption spectra + ( + [4, 1, 1], + [], + np.array([np.eye(3) * k_point * (1 + 1j) for k_point in range(1, 5)]), + None, + ), + # Even if we define `Frequency`, we cannot extract absorption spectra if `value` depends on `KMesh` + ( + [4, 1, 1], + [0, 1, 2, 3, 4], + np.array( + [ + [ + np.eye(3) * k_point * (1 + 1j) + + np.eye(3) * freq_point * 0.5j + for freq_point in range(5) + ] + for k_point in range(1, 5) + ] + ), + None, + ), + # Valid case: `value` does not depend on `KMesh` and we can extract absorption spectra + ( + [], + [0, 1, 2, 3, 4], + np.array([np.eye(3) * freq_point * 0.5j for freq_point in range(5)]), + [0.0, 0.5, 1.0, 1.5, 2.0], + ), + ], + ) + def test_extract_absorption_spectra( + self, + kmesh_grid: Optional[List[int]], + frequencies_points: Optional[List[float]], + value: Optional[np.ndarray], + result: Optional[List[float]], + ): + """ + Test the `extract_absorption_spectra` method. The `result` in the last valid case corresponds to the imaginary part of + the diagonal of the `Permittivity.value` for each frequency point. + """ + permittivity = Permittivity() + + # Generating `KMesh` numerical settings section + simulation = generate_k_space_simulation(grid=kmesh_grid) + k_mesh_settings = simulation.model_method[0].numerical_settings[0].k_mesh[0] + k_mesh_settings.label = 'q-mesh' + k_mesh_settings.center = 'Monkhorst-Pack' + k_mesh_settings.points, _ = k_mesh_settings.resolve_points_and_offset(logger) + if kmesh_grid is not None and len(kmesh_grid) > 0: + kmesh_variables = KMesh(points=k_mesh_settings) + permittivity.variables.append(kmesh_variables) + + # Adding `Frequency` if defined + if frequencies_points is not None and len(frequencies_points) > 0: + frequencies = Frequency(points=frequencies_points) + permittivity.variables.append(frequencies) + + if permittivity.variables is not None and len(permittivity.variables) > 0: + permittivity.value = value + + absorption_spectra = permittivity.extract_absorption_spectra(logger) + if absorption_spectra is not None: + assert len(absorption_spectra) == 3 + spectrum = absorption_spectra[1] + assert spectrum.rank == [] + assert spectrum.axis == 'yy' + assert len(spectrum.value) == len(permittivity.variables[0].points) + assert np.allclose(spectrum.value, result) + else: + assert absorption_spectra == result diff --git a/tests/test_spectral_profile.py b/tests/test_spectral_profile.py index 4316ecf7..14a55034 100644 --- a/tests/test_spectral_profile.py +++ b/tests/test_spectral_profile.py @@ -32,9 +32,10 @@ from nomad_simulations.properties import ( SpectralProfile, ElectronicDensityOfStates, - XASSpectra, + AbsorptionSpectrum, + XASSpectrum, ) -from nomad_simulations.variables import Temperature, Energy2 as Energy +from nomad_simulations.variables import Energy2 as Energy class TestSpectralProfile: @@ -73,19 +74,6 @@ def test_default_quantities(self): assert electronic_dos.name == 'ElectronicDensityOfStates' assert electronic_dos.rank == [] - def test_get_energy_points(self): - """ - Test the `_get_energy_points` private method. We test here that the function indeed returns the points from a `Energy` variable. - """ - electronic_dos = ElectronicDensityOfStates() - electronic_dos.variables = [ - Temperature(points=list(range(-3, 4)) * ureg.kelvin) - ] - assert electronic_dos._get_energy_points(logger) is None - electronic_dos.variables.append(Energy(points=list(range(-3, 4)) * ureg.joule)) - energies = electronic_dos._get_energy_points(logger) - assert (energies.magnitude == np.array([-3, -2, -1, 0, 1, 2, 3])).all() - def test_resolve_energies_origin(self): """ Test the `resolve_energies_origin` method. @@ -249,20 +237,36 @@ def test_normalize(self): pass -class TestXASSpectra: +class TestAbsorptionSpectrum: + """ + Test the `AbsorptionSpectrum` class defined in `properties/spectral_profile.py`. + """ + + # ! Include this initial `test_default_quantities` method when testing your PhysicalProperty classes + def test_default_quantities(self): + """ + Test the default quantities assigned when creating an instance of the `AbsorptionSpectrum` class. + """ + absorption_spectrum = AbsorptionSpectrum() + assert absorption_spectrum.iri is None # Add iri when available + assert absorption_spectrum.name == 'AbsorptionSpectrum' + assert absorption_spectrum.rank == [] + + +class TestXASSpectrum: """ - Test the `XASSpectra` class defined in `properties/spectral_profile.py`. + Test the `XASSpectrum` class defined in `properties/spectral_profile.py`. """ # ! Include this initial `test_default_quantities` method when testing your PhysicalProperty classes def test_default_quantities(self): """ - Test the default quantities assigned when creating an instance of the `XASSpectra` class. + Test the default quantities assigned when creating an instance of the `XASSpectrum` class. """ - xas_spectra = XASSpectra() - assert xas_spectra.iri is None # Add iri when available - assert xas_spectra.name == 'XASSpectra' - assert xas_spectra.rank == [] + xas_spectrum = XASSpectrum() + assert xas_spectrum.iri is None # Add iri when available + assert xas_spectrum.name == 'XASSpectrum' + assert xas_spectrum.rank == [] @pytest.mark.parametrize( 'xanes_energies, exafs_energies, xas_values', @@ -284,19 +288,19 @@ def test_generate_from_contributions( """ Test the `generate_from_contributions` method. """ - xas_spectra = XASSpectra() + xas_spectrum = XASSpectrum() if xanes_energies is not None: - xanes_spectra = SpectralProfile() - xanes_spectra.variables = [Energy(points=xanes_energies * ureg.joule)] - xanes_spectra.value = [0.5, 0.1, 0.3] - xas_spectra.xanes_spectra = xanes_spectra + xanes_spectrum = AbsorptionSpectrum() + xanes_spectrum.variables = [Energy(points=xanes_energies * ureg.joule)] + xanes_spectrum.value = [0.5, 0.1, 0.3] + xas_spectrum.xanes_spectrum = xanes_spectrum if exafs_energies is not None: - exafs_spectra = SpectralProfile() - exafs_spectra.variables = [Energy(points=exafs_energies * ureg.joule)] - exafs_spectra.value = [0.2, 0.4, 0.6] - xas_spectra.exafs_spectra = exafs_spectra - xas_spectra.generate_from_contributions(logger) - if xas_spectra.value is not None: - assert (xas_spectra.value == xas_values).all() + exafs_spectrum = AbsorptionSpectrum() + exafs_spectrum.variables = [Energy(points=exafs_energies * ureg.joule)] + exafs_spectrum.value = [0.2, 0.4, 0.6] + xas_spectrum.exafs_spectrum = exafs_spectrum + xas_spectrum.generate_from_contributions(logger) + if xas_spectrum.value is not None: + assert (xas_spectrum.value == xas_values).all() else: - assert xas_spectra.value == xas_values + assert xas_spectrum.value == xas_values diff --git a/tests/test_utils.py b/tests/test_utils.py index bdbe6e1c..5a498952 100644 --- a/tests/test_utils.py +++ b/tests/test_utils.py @@ -16,13 +16,20 @@ # limitations under the License. # -from . import logger +import pytest +import numpy as np + +from nomad.units import ureg from nomad_simulations.utils import ( get_sibling_section, is_not_representative, + get_variables, ) from nomad_simulations.model_system import ModelSystem, AtomicCell, Symmetry +from nomad_simulations.variables import Temperature, Energy2 as Energy + +from . import logger def test_get_sibling_section(): @@ -59,3 +66,28 @@ def test_is_not_representative(): # ! Missing test for RusselSandersState (but this class will probably be deprecated) + + +@pytest.mark.parametrize( + 'variables, result, result_length', + [ + (None, [], 0), + ([], [], 0), + ([Temperature()], [], 0), + ([Temperature(), Energy(n_points=4)], [Energy(n_points=4)], 1), + ( + [Temperature(), Energy(n_points=2), Energy(n_points=10)], + [Energy(n_points=2), Energy(n_points=10)], + 2, + ), + # TODO add testing when we have variables which inherit from another variable + ], +) +def test_get_variables(variables: list, result: list, result_length: int): + """ + Test the `get_variables` utility function + """ + energies = get_variables(variables, Energy) + assert len(energies) == result_length + for i, energy in enumerate(energies): # asserting energies == result does not work + assert energy.n_points == result[i].n_points