import logging
from abc import ABC, abstractmethod
from typing import Union
import numpy as np
from ..types import Array, ComplexArray, Path_str
from ..units import au_to_eV, eV_to_au
from ..utilities import ClassFormatter
logger = logging.getLogger(__name__)
[docs]
class DielectricFunction(ABC):
"""Objects of this class hold different representations of dielectric functions.
Parameters
----------
name
Name of dielectric function.
"""
def __init__(self,
name: str = 'DielectricFunction') -> None:
logger.debug(f'Entering {self.__class__.__name__}.__init__')
self._name = name
[docs]
def __call__(self, frequencies: Array) -> np.ndarray:
"""Returns the dielectric function at the given frequencies
Parameters
----------
frequencies
List of frequencies at which the dielectric function is returned;
in units of eV.
"""
freq_w = np.array(frequencies, dtype=float)
freq_w *= eV_to_au
return self._eval_at(freq_w)
@abstractmethod
def _eval_at(self, freq_w: np.ndarray) -> np.ndarray:
"""Returns the dielectric function at the given frequencies
Parameters
----------
frequencies
List of frequencies at which the dielectric function is returned;
in atomic units.
"""
raise NotImplementedError()
@abstractmethod
def __str__(self) -> str:
raise NotImplementedError()
@property
def name(self) -> str:
"""Name of dielectric function."""
return self._name
class AnalyticDielectricFunction(DielectricFunction):
"""Objects of this class hold different representations of dielectric functions that
have analytic form. The dielectric function can thus be evaluated at any frequency
without interpolation.
Parameters
----------
name
Name of dielectric function.
"""
def __init__(self,
name: str = 'AnalyticDielectricFunction') -> None:
logger.debug(f'Entering {self.__class__.__name__}.__init__')
super().__init__(name=name)
[docs]
class NumericDielectricFunction(DielectricFunction):
"""Objects of this class hold different representations of dielectric functions sampled
on a discrete grid of frequencies. The dielectric function can be evaluated at any frequency
within the supplied grid by linear interpolation.
Parameters
----------
frequencies
List of frequencies at which the dielectric function :attr:`dielectric_function`
is sampled; units are eV by default, optionally atomic units
(see :attr:`units`).
The list of frequencies and the corresponding dielectric function are sorted if the former
is not monotonically increasing.
dielectric_function
Complex dielectric function sampled at the frequencies :attr:`frequencies`.
units
`eV` to specify :attr:`freq_w` in eV, or `au` to specify in atomic units.
This parameter determines whether conversion should be performed during
initialization and has no effect on instance methods and variables.
name
Name of dielectric function.
Raises
------
ValueError
If :attr:`frequencies` or :attr:`dielectric_function` have zero length.
ValueError
If there is a mismatch in the length of :attr:`frequencies` and :attr:`dielectric_function`.
ValueError
If :attr:`units` is not one of the allowed values `eV` or `au`.
Examples
--------
The following snippet illustrates how a NumericDielectricFunction object can be constructed.
Here, the analytic Drude dielectric function of the free-electron gas is used.
To construct a dielectric function object from a numeric dielectric function the function
:func:`~strongcoca.response.read_dielec_function_file` can be used.
>>> import numpy as np
>>> from strongcoca.response.dielectric import NumericDielectricFunction
>>>
>>> omega_p = 3
>>> gamma = 0.5
>>> freq_w = np.linspace(0.1, 10)
>>> eps_w = 1 - omega_p**2 / (freq_w**2 + 1j*gamma*freq_w)
>>> dielec = NumericDielectricFunction(freq_w, eps_w, name='Drude-function')
>>> print(dielec)
---------------- NumericDielectricFunction -----------------
name : Drude-function
Frequency grid : 50 points between 0.10 and 10.00eV
frequencies (eV) : [ 0.1 0.30204 0.50408 ...
9.59592 9.79796 10. ]
dielectric_function : [-33.61538+1.73077e+02j
-25.37528+4.36618e+01j
-16.85366+1.77091e+01j ...
0.90253+5.07897e-03j
0.90649+4.77173e-03j
0.91022+4.48878e-03j]
>>> dense_freq_w = np.linspace(3, 4)
>>> print(dielec(dense_freq_w))
[0.02425453+0.16310249j 0.0367996 +0.15996402j 0.04934466+0.15682555j
0.06188972+0.15368708j 0.07443478+0.15054861j 0.08697984+0.14741014j
...
"""
def __init__(self,
frequencies: Union[Array, float],
dielectric_function: Union[ComplexArray, complex],
units: str = 'eV',
name: str = 'NumericDielectricFunction') -> None:
logger.debug(f'Entering {self.__class__.__name__}.__init__')
super().__init__(name=name)
if isinstance(frequencies, float):
frequencies = [frequencies]
if isinstance(dielectric_function, complex) or isinstance(dielectric_function, float):
dielectric_function = [dielectric_function]
freq_w = np.array(frequencies)
eps_w = np.array(dielectric_function, dtype=complex)
if units == 'eV':
freq_w = freq_w * eV_to_au
elif units != 'au':
raise ValueError(f"units has to be 'eV' or 'au', not '{units}'")
if len(freq_w) != len(eps_w):
raise ValueError(f'length of frequencies {len(freq_w)} must match length of '
f'dielectric function {len(eps_w)}')
if len(freq_w) == 0:
raise ValueError('length of frequencies and dielectric function must not be zero')
# np.diff works if input has at least length 2
if len(freq_w) > 1 and not np.all(np.diff(freq_w) > 0):
sort_w = np.argsort(freq_w)
freq_w = freq_w[sort_w]
eps_w = eps_w[sort_w]
self._freq_w = freq_w
self._eps_w = eps_w
def _eval_at(self, freq_w: np.ndarray) -> np.ndarray:
"""Returns the dielectric function at the given frequencies by interpolating the numeric
dielectric function.
Parameters
----------
frequencies
List of frequencies at which the dielectric function is returned;
in units of eV.
Raises
------
ValueError
If any component of :attr:`frequencies` is outside of the interval of the numeric
dielectric function.
"""
if np.min(freq_w) < self._freq_w[0]:
raise ValueError('Minimum of array frequencies is outside the lower bound of the'
'numeric dielectric function')
if np.max(freq_w) > self._freq_w[-1]:
raise ValueError('Maximum of array frequencies is outside the upper bound of the'
'numeric dielectric function')
eps_w = np.interp(freq_w, self._freq_w, self._eps_w)
return eps_w # type: ignore
def __str__(self) -> str:
fmt = ClassFormatter(self, pad=20)
fmt.append_class_name()
fmt.append_attr('name')
fmt.append_attr('Frequency grid', f'{self.n} points between {self.minfreq:.2f} and '
f'{self.maxfreq:.2f}eV')
fmt.append_attr('frequencies', unit='eV')
fmt.append_attr('dielectric_function')
return fmt.to_string()
@property
def n(self) -> int:
"""number of samples for the dielectric function
:attr:`~strongcoca.response.NumericDielectricFunction.dielectric_function`."""
return len(self._eps_w)
@property
def minfreq(self) -> float:
"""Lowest frequency :attr:`~strongcoca.response.NumericDielectricFunction.frequencies`
for which there is a value of the dielectric function; in units of eV."""
m = self._freq_w[0] * au_to_eV
assert isinstance(m, float)
return m
@property
def maxfreq(self) -> float:
"""Highest frequency :attr:`~strongcoca.response.NumericDielectricFunction.frequencies`
for which there is a value of the dielectric function; in units of eV."""
m = self._freq_w[-1] * au_to_eV
assert isinstance(m, float)
return m
@property
def frequencies(self) -> np.ndarray:
"""Frequencies corresponding to dielectric function
:attr:`~strongcoca.response.NumericDielectricFunction.dielectric_function` in units of eV;
array with shape {n}."""
return self._freq_w * au_to_eV # type: ignore
@property
def dielectric_function(self) -> np.ndarray:
"""Complex unitless dielectric function at frequencies
:attr:`~strongcoca.response.MieGansResponse.frequencies`; array with shape {n}."""
return self._eps_w
[docs]
class DrudeDielectricFunction(AnalyticDielectricFunction):
r"""Objects of this class represent the Drude dielectric function, which is defined for the
entire positive frequency domain. The function is implemented as
.. math::
\varepsilon_r(\omega) = 1 - \frac{\omega_p^2}{\omega^2 + i\gamma\omega}
where the plasma frequency :math:`\omega_p` determines where the real part of the dielectric
function is zero and :math:`\gamma` is a broadening.
Parameters
----------
plasma_frequency
The plasma frequency :math:`\omega_p`; units are eV by default, optionally atomic units
(see :attr:`units`).
broadening
Broadening of the dielectric function :math:`\gamma` units are eV by default, optionally
atomic units (see :attr:`units`).
name
Name of the dielectric function.
units
`eV` to specify :attr:`plasma_frequency` and :attr:`broadening` in eV, or `au` to specify
in atomic units.
This parameter determines whether conversion should be performed during
initialization and has no effect on instance methods and variables.
"""
def __init__(self,
plasma_frequency: float,
broadening: float,
units: str = 'eV',
name: str = 'DrudeDielectricFunction') -> None:
logger.debug(f'Entering {self.__class__.__name__}.__init__')
super().__init__(name=name)
omega_pl = plasma_frequency
gamma = broadening
if units == 'eV':
omega_pl *= eV_to_au
gamma *= eV_to_au
elif units != 'au':
raise ValueError(f"units has to be 'eV' or 'au', not '{units}'")
self._omega_pl = omega_pl
self._gamma = gamma
def _eval_at(self, freq_w: np.ndarray) -> np.ndarray:
eps_w = 1 - self._omega_pl**2 / (freq_w**2 + 1j*self._gamma*freq_w)
return eps_w # type: ignore
def __str__(self) -> str:
fmt = ClassFormatter(self, pad=15)
fmt.append_class_name()
fmt.append_attr('name')
fmt.append_attr('plasma_frequency', unit='eV')
fmt.append_attr('broadening', unit='eV')
return fmt.to_string()
@property
def plasma_frequency(self) -> float:
"""plasma frequency."""
return self._omega_pl * au_to_eV
@property
def broadening(self) -> float:
"""Broadening parameter."""
return self._gamma * au_to_eV
[docs]
def read_dielec_function_file(dielectric_function_file: Path_str,
frequency_column: int = 0,
dielectric_function_real_column: int = 1,
dielectric_function_imag_column: int = 2,
units: str = 'eV',
name: str = 'NumericDielectricFunction from file',
**np_loadtxt_kwargs,
) -> NumericDielectricFunction:
"""Parses a textfile with columns describing the dielectric function of a material.
Parameters
----------
dielectric_function_file
Path to the file containing the dielectic function.
frequency_column
Index of column containing the frequency vector in units of eV.
dielectric_function_real_column
Index of column containing the real part of the unitless dielectric function.
dielectric_function_imag_column
Index of column containing the imaginary part of the unitless dielectric function.
units
Units of the frequency column;
passed to :class:`strongcoca.response.NumericDielectricFunction`.
name
Passed to :class:`~strongcoca.response.NumericDielectricFunction`.
np_loadtxt_kwargs
Parameters passed to :func:`numpy.loadtxt`; useful parameters are `delimiter` or `comments`.
"""
freq_w, df_real_w, df_imag_w = np.loadtxt(dielectric_function_file,
usecols=[frequency_column,
dielectric_function_real_column,
dielectric_function_imag_column],
unpack=True,
**np_loadtxt_kwargs)
df_w = df_real_w + 1j*df_imag_w
df = NumericDielectricFunction(freq_w, df_w, units=units, name=name)
return df
[docs]
def read_dielec_function_from_refractive_index_file(
dielectric_function_file: Path_str,
frequency_column: int = 0,
refractive_index_column: int = 1,
extinction_coefficient_column: int = 2,
units: str = 'eV',
name: str = 'NumericDielectricFunction from file',
**np_loadtxt_kwargs,
) -> NumericDielectricFunction:
r"""Parses a textfile with columns describing the refractive index and extinction coefficient
of a material.
The refractive index :math:`n(\omega)` and extinction coefficient :math:`k(\omega)` are
related to the complex dielectric function :math:`\varepsilon_r(\omega)` by
.. math::
\varepsilon_r(\omega) = n^2(\omega) - k^2(\omega) + 2jn(\omega)k(\omega)
Parameters
----------
dielectric_function_file
Path to the file containing the dielectic function.
frequency_column
Index of column containing the frequency vector in units of eV.
refractive_index_column
Index of column containing the unitless refractive index.
extinction_coefficient_column
Index of column containing the unitless extinction coefficient.
units
Units of the frequency column;
passed to :class:`~strongcoca.response.NumericDielectricFunction`.
name
Passed to :class:`~strongcoca.response.NumericDielectricFunction`.
np_loadtxt_kwargs
Parameters passed to :func:`numpy.loadtxt`; useful parameters are `delimiter` or `comments`.
"""
freq_w, n_w, k_w = np.loadtxt(dielectric_function_file,
usecols=[frequency_column,
refractive_index_column,
extinction_coefficient_column],
unpack=True,
**np_loadtxt_kwargs)
df_w = n_w**2 - k_w**2 + 2j*n_w*k_w
df = NumericDielectricFunction(freq_w, df_w, units=units, name=name)
return df
