Squashing stuff to avoid issues when rebasing...

src/nomad_simulations/numerical_settings.py

@@ -462,6 +462,99 @@ def resolve_k_line_density(
                 return k_line_density
         return None
 
+    def resolve_high_symmetry_points(
+        self,
+        model_systems: List[ModelSystem],
+        logger: BoundLogger,
+        eps: float = 3e-3,
+    ) -> Optional[dict]:
+        """
+        Resolves the `high_symmetry_points` of the `KMesh` from the list of `ModelSystem`. This method
+        relies on using the `ModelSystem` information in the sub-sections `Symmetry` and `AtomicCell`, and uses
+        the ASE package to extract the special (high symmetry) points information.
+
+        Args:
+            model_systems (List[ModelSystem]): The list of `ModelSystem` sections.
+            logger (BoundLogger): The logger to log messages.
+            eps (float, optional): Tolerance factor to define the `lattice` ASE object. Defaults to 3e-3.
+
+        Returns:
+            (Optional[dict]): The resolved `high_symmetry_points` of the `KMesh`.
+        """
+        # Extracting `bravais_lattice` from `ModelSystem.symmetry` section and `ASE.cell` from `ModelSystem.cell`
+        lattice = None
+        for model_system in model_systems:
+            # General checks to proceed with normalization
+            if is_not_representative(model_system, logger):
+                continue
+            if model_system.symmetry is None:
+                logger.warning('Could not find `ModelSystem.symmetry`.')
+                continue
+            bravais_lattice = [symm.bravais_lattice for symm in model_system.symmetry]
+            if len(bravais_lattice) != 1:
+                logger.warning(
+                    'Could not uniquely determine `bravais_lattice` from `ModelSystem.symmetry`.'
+                )
+                continue
+            bravais_lattice = bravais_lattice[0]
+
+            if model_system.cell is None:
+                logger.warning('Could not find `ModelSystem.cell`.')
+                continue
+            prim_atomic_cell = None
+            for atomic_cell in model_system.cell:
+                if atomic_cell.type == 'primitive':
+                    prim_atomic_cell = atomic_cell
+                    break
+            if prim_atomic_cell is None:
+                logger.warning(
+                    'Could not find the primitive `AtomicCell` under `ModelSystem.cell`.'
+                )
+                continue
+            # function defined in AtomicCell
+            atoms = prim_atomic_cell.to_ase_atoms(logger)
+            cell = atoms.get_cell()
+            lattice = cell.get_bravais_lattice(eps)
+            break  # only cover the first representative `ModelSystem`
+
+        # Checking if `bravais_lattice` and `lattice` are defined
+        if lattice is None:
+            logger.warning(
+                'Could not resolve `bravais_lattice` and `lattice` ASE object from the `ModelSystem`.'
+            )
+            return None
+
+        # Non-conventional ordering testing for certain lattices:
+        if bravais_lattice in ['oP', 'oF', 'oI', 'oS']:
+            a, b, c = lattice.a, lattice.b, lattice.c
+            assert a < b
+            if bravais_lattice != 'oS':
+                assert b < c
+        elif bravais_lattice in ['mP', 'mS']:
+            a, b, c = lattice.a, lattice.b, lattice.c
+            alpha = lattice.alpha * np.pi / 180
+            assert a <= c and b <= c  # ordering of the conventional lattice
+            assert alpha < np.pi / 2
+
+        # Extracting the `high_symmetry_points` from the `lattice` object
+        special_points = lattice.get_special_points()
+        if special_points is None:
+            logger.warning(
+                'Could not find `lattice.get_special_points()` from the ASE package.'
+            )
+            return None
+        high_symmetry_points = {}
+        for key, value in lattice.get_special_points().items():
+            if key == 'G':
+                key = 'Gamma'
+            if bravais_lattice == 'tI':
+                if key == 'S':
+                    key = 'Sigma'
+                elif key == 'S1':
+                    key = 'Sigma1'
+            high_symmetry_points[key] = list(value)
+        return high_symmetry_points
+
     def normalize(self, archive, logger) -> None:
         super().normalize(archive, logger)
 
@@ -679,6 +772,7 @@ def resolve_points(
                 'The `reciprocal_lattice_vectors` are not passed as an input.'
             )
             return None
+
         # Check if `points_norm` is a list and convert it to a numpy array
         if isinstance(points_norm, list):
             points_norm = np.array(points_norm)

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)
 
     # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # #
     # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # # #

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

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