diff --git a/docs/coverage-badge.svg b/docs/coverage-badge.svg index 8691f97..83b2e82 100644 --- a/docs/coverage-badge.svg +++ b/docs/coverage-badge.svg @@ -1 +1 @@ -coverage: 93.46%coverage93.46% +coverage: 93.22%coverage93.22% diff --git a/docs/source/generated/ramannoodle.polarizability.rst b/docs/source/generated/ramannoodle.polarizability.rst index 555192f..d27d7ab 100644 --- a/docs/source/generated/ramannoodle.polarizability.rst +++ b/docs/source/generated/ramannoodle.polarizability.rst @@ -4,6 +4,14 @@ ramannoodle.polarizability package Submodules ---------- +ramannoodle.polarizability.art module +------------------------------------- + +.. automodule:: ramannoodle.polarizability.art + :members: + :undoc-members: + :show-inheritance: + ramannoodle.polarizability.interpolation module ----------------------------------------------- diff --git a/docs/tests-badge.svg b/docs/tests-badge.svg index 9611b98..5ee71a0 100644 --- a/docs/tests-badge.svg +++ b/docs/tests-badge.svg @@ -1 +1 @@ -tests: 87tests87 +tests: 104tests104 diff --git a/ramannoodle/polarizability/art.py b/ramannoodle/polarizability/art.py new file mode 100644 index 0000000..e8c5a11 --- /dev/null +++ b/ramannoodle/polarizability/art.py @@ -0,0 +1,191 @@ +"""Polarizability model based on atomic Raman tensors (ARTs).""" + +from pathlib import Path + +import numpy as np +from numpy.typing import NDArray + +from .interpolation import InterpolationModel +from ..exceptions import ( + get_type_error, + get_shape_error, + verify_ndarray_shape, + InvalidDOFException, +) + + +class ARTModel(InterpolationModel): + """Polarizability model based on atomic Raman tensors (ARTs). + + This model is a flavor of InterpolationModel with added restrictions. The degrees + of freedom (DOF) in ARTModel are displacements of individual atoms. 1D + interpolations are used, with only two amplitudes allowed per DOF. Atomic + displacements can be in any direction, so long as all three directions for a given + atom are orthogonal. + + All assumptions from InterpolationModel apply. + + .. note:: + Degrees of freedom cannot (and should not) be added using `add_dof` and + `add_dof_from_files`. Instead, use `add_art` and `add_art_from_files`. + + .. warning:: + Due to the way InterpolationModel works, an ARTModel's predicted polarizability + tensor of the equilibrium structure may be offset slightly. This is a fixed + offset and therefore has no influence on Raman spectra calculations. + + Parameters + ---------- + structural_symmetry + equilibrium_polarizability + 2D array with shape (3,3) giving polarizability of system at equilibrium. This + would usually correspond to the minimum energy structure. + + """ + + def add_dof( # pylint: disable=too-many-arguments + self, + displacement: NDArray[np.float64], + amplitudes: NDArray[np.float64], + polarizabilities: NDArray[np.float64], + interpolation_order: int, + include_equilibrium_polarizability: bool = True, + ) -> None: + """Disable add_dof. + + :meta private: + """ + raise AttributeError("'ARTModel' object has no attribute 'add_dof'") + + def add_dof_from_files( + self, + filepaths: str | Path | list[str] | list[Path], + file_format: str, + interpolation_order: int, + ) -> None: + """Disable add_dof_from_files. + + :meta private: + """ + raise AttributeError("'ARTModel' object has no attribute 'add_dof_from_files'") + + def add_art( + self, + atom_index: int, + direction: NDArray[np.float64], + amplitudes: NDArray[np.float64], + polarizabilities: NDArray[np.float64], + ) -> None: + """Add an atomic Raman tensor (ART). + + Specification of an ART requires an atom index, displacement direction, + displacement amplitudes, and corresponding known polarizabilities for each + amplitude. Alongside ART's related to the specified ART by symmetry are added + automatically. + + Parameters + ---------- + atom_index + Index of atom, consistent with the StructuralSymmetry used to initialize + the model. + direction + 1D array with shape (3,). Must be orthogonal to all previously added ARTs + belonging to specified atom. + amplitudes + 1D array of length 1 or 2 containing amplitudes in angstroms. Duplicate + amplitudes are not allowed, including symmetrically equivalent + amplitudes. + polarizabilities + 3D array with shape (1,3,3) or (2,3,3) containing known polarizabilities for + each amplitude. + + Raises + ------ + InvalidDOFException + Provided ART was invalid. + + """ + if not isinstance(atom_index, int): + raise get_type_error("atom_index", atom_index, "int") + try: + if amplitudes.shape not in [(1,), (2,)]: + raise get_shape_error("amplitudes", amplitudes, "(1,) or (2,)") + except AttributeError as exc: + raise get_type_error("amplitudes", amplitudes, "ndarray") from exc + verify_ndarray_shape( + "polarizabilities", polarizabilities, (amplitudes.size, 3, 3) + ) + + displacement = self._structural_symmetry.get_fractional_positions() * 0 + try: + displacement[atom_index] = direction / np.linalg.norm(direction * 10.0) + except TypeError as exc: + raise get_type_error("direction", direction, "ndarray") from exc + except ValueError as exc: + raise get_shape_error("direction", direction, "(3,)") from exc + + super().add_dof( + displacement, + amplitudes, + polarizabilities, + 1, + include_equilibrium_polarizability=False, + ) + + def add_art_from_files( + self, + filepaths: str | Path | list[str] | list[Path], + file_format: str, + ) -> None: + """Add a atomic Raman tensor (ART) from file(s). + + Required directions, amplitudes, and polarizabilities are automatically + determined from provided files. Files should be chosen wisely such that the + resulting ARTs are valid under the restrictions set by `add_art`. + + Parameters + ---------- + filepaths + file_format + supports: "outcar" + + Raises + ------ + FileNotFoundError + File could not be found. + InvalidDOFException + ART assembled from supplied files was invalid (see get_art) + + """ + displacements, amplitudes, polarizabilities = super()._read_dof( + filepaths, file_format + ) + # Check whether only one atom is displaced. + _displacement = np.copy(displacements[0]) + atom_index = int(np.argmax(np.sum(_displacement**2, axis=1))) + _displacement[atom_index] = np.zeros(3) + if not np.isclose(_displacement, 0.0, atol=1e-6).all(): + raise InvalidDOFException("multiple atoms displaced simultaneously") + + basis_vectors_to_add, interpolation_xs, interpolation_ys = super()._get_dof( + displacements[0], amplitudes, polarizabilities, False + ) + + num_amplitudes = len(interpolation_xs[0]) + if num_amplitudes != 2: + raise InvalidDOFException( + f"wrong number of amplitudes: {num_amplitudes} != 2" + ) + + super()._construct_and_add_interpolations( + basis_vectors_to_add, interpolation_xs, interpolation_ys, 1 + ) + + def get_status_dict(self) -> dict[str, int]: + """Return information on model status. + + Relevant information includes which ARTs are required for each atom and whether + or not these ARTs have been specified. + + """ + raise NotImplementedError diff --git a/test/data/TiO2/known_art_spectrum.npz b/test/data/TiO2/known_art_spectrum.npz new file mode 100644 index 0000000..0fc48c0 Binary files /dev/null and b/test/data/TiO2/known_art_spectrum.npz differ diff --git a/test/tests/test_art.py b/test/tests/test_art.py new file mode 100644 index 0000000..d749ff4 --- /dev/null +++ b/test/tests/test_art.py @@ -0,0 +1,211 @@ +"""Testing for ARTModel.""" + +from typing import Type + +import numpy as np +from numpy.typing import NDArray + +import pytest + +from ramannoodle.polarizability.art import ARTModel +from ramannoodle.exceptions import InvalidDOFException +from ramannoodle.symmetry import StructuralSymmetry + +# pylint: disable=protected-access +# pylint: disable=too-many-arguments + + +@pytest.mark.parametrize( + "outcar_symmetry_fixture,atom_index, direction, amplitudes, known_dof_added", + [ + ( + "test/data/STO_RATTLED_OUTCAR", + 0, + np.array([1, 0, 0]), + np.array([-0.01, 0.01]), + 1, + ), + ("test/data/TiO2/phonons_OUTCAR", 0, np.array([1, 0, 0]), np.array([0.01]), 72), + ], + indirect=["outcar_symmetry_fixture"], +) +def test_add_art( + outcar_symmetry_fixture: StructuralSymmetry, + atom_index: int, + direction: NDArray[np.float64], + amplitudes: NDArray[np.float64], + known_dof_added: int, +) -> None: + """Test add_art (normal).""" + symmetry = outcar_symmetry_fixture + model = ARTModel(symmetry, np.zeros((3, 3))) + model.add_art(atom_index, direction, amplitudes, np.zeros((amplitudes.size, 3, 3))) + assert len(model._cartesian_basis_vectors) == known_dof_added + assert np.isclose(np.linalg.norm(model._cartesian_basis_vectors[0]), 1) + + +@pytest.mark.parametrize( + "outcar_symmetry_fixture,atom_indexes,directions,amplitudes,polarizabilities," + "exception_type,in_reason", + [ + ( + "test/data/STO_RATTLED_OUTCAR", + [0], + np.array([[1, 0, 0]]), + np.array([0.1]), + np.zeros((1, 3, 3)), + InvalidDOFException, + "insufficient points", + ), + ( + "test/data/STO_RATTLED_OUTCAR", + [0], + np.array([[1, 0, 0]]), + np.array([-0.1, 0.1]), + np.zeros((2, 4, 3)), + ValueError, + "polarizabilities has wrong shape", + ), + ( + "test/data/STO_RATTLED_OUTCAR", + [0], + np.array([[1, 0, 0]]), + np.array([-0.1, 0.1]), + np.zeros((4, 3, 3)), + ValueError, + "polarizabilities has wrong shape", + ), + ( + "test/data/STO_RATTLED_OUTCAR", + [0, 0], + np.array([[1, 0, 0], [1, 0, 0]]), + np.array([-0.1, 0.1]), + np.zeros((2, 3, 3)), + InvalidDOFException, + "not orthogonal", + ), + ( + "test/data/STO_RATTLED_OUTCAR", + [0, 0], + np.array([[1, 0, 0], [1, 0, 0]]), + [-0.1, 0.1], + np.zeros((2, 3, 3)), + TypeError, + "amplitudes should have type ndarray, not list", + ), + ( + "test/data/STO_RATTLED_OUTCAR", + [0, 0], + np.array([[1, 0, 0], [1, 0, 0]]), + np.array([-0.1, 0.1]), + list(np.zeros((2, 3, 3))), + TypeError, + "polarizabilities should have type ndarray, not list", + ), + ( + "test/data/STO_RATTLED_OUTCAR", + [[0, 0]], + np.array([[1, 0, 0], [1, 0, 0]]), + np.array([-0.1, 0.1]), + np.zeros((2, 3, 3)), + TypeError, + "atom_index should have type int, not list", + ), + ( + "test/data/STO_RATTLED_OUTCAR", + [0], + np.array([[1, 0, 0], [1, 0, 0]]), + np.array([0.1, 0.1]), + np.zeros((2, 3, 3)), + InvalidDOFException, + "due to symmetry, amplitude 0.1 should not be specified", + ), + ], + indirect=["outcar_symmetry_fixture"], +) +def test_add_art_exception( + outcar_symmetry_fixture: StructuralSymmetry, + atom_indexes: list[int], + directions: NDArray[np.float64], + amplitudes: NDArray[np.float64], + polarizabilities: NDArray[np.float64], + exception_type: Type[Exception], + in_reason: str, +) -> None: + """Test add_art (exception).""" + symmetry = outcar_symmetry_fixture + model = ARTModel(symmetry, np.zeros((3, 3))) + with pytest.raises(exception_type) as error: + for atom_index, direction in zip(atom_indexes, directions): + model.add_art(atom_index, direction, amplitudes, polarizabilities) + + assert in_reason in str(error.value) + + +@pytest.mark.parametrize( + "outcar_symmetry_fixture,outcar_files,exception_type,in_reason", + [ + ( + "test/data/STO_RATTLED_OUTCAR", + ["test/data/TiO2/Ti5_0.1x_eps_OUTCAR"], + InvalidDOFException, + "incompatible outcar", + ), + ( + "test/data/TiO2/phonons_OUTCAR", + [ + "test/data/TiO2/Ti5_0.1x_eps_OUTCAR", + "test/data/TiO2/Ti5_0.2x_eps_OUTCAR", + ], + InvalidDOFException, + "wrong number of amplitudes: 4 != 2", + ), + ( + "test/data/TiO2/phonons_OUTCAR", + [ + "test/data/TiO2/Ti5_0.1x_eps_OUTCAR", + "test/data/TiO2/Ti5_m0.1x_eps_OUTCAR", + ], + InvalidDOFException, + "wrong number of amplitudes: 4 != 2", + ), + ( + "test/data/TiO2/phonons_OUTCAR", + [ + "this_outcar_does_not_exist", + ], + FileNotFoundError, + "No such file or directory", + ), + ( + "test/data/TiO2/phonons_OUTCAR", + [ + "test/data/TiO2/Ti5_0.1x_eps_OUTCAR", + "test/data/TiO2/Ti5_0.1y_eps_OUTCAR", + ], + InvalidDOFException, + "is not collinear", + ), + ( + "test/data/TiO2/phonons_OUTCAR", + [ + "test/data/TiO2/O43_0.1z_eps_OUTCAR", + ], + InvalidDOFException, + "wrong number of amplitudes: 1 != 2", + ), + ], + indirect=["outcar_symmetry_fixture"], +) +def test_add_art_from_files_exception( + outcar_symmetry_fixture: StructuralSymmetry, + outcar_files: list[str], + exception_type: Type[Exception], + in_reason: str, +) -> None: + """Test add_dof_from_files (exception).""" + symmetry = outcar_symmetry_fixture + model = ARTModel(symmetry, np.zeros((3, 3))) + with pytest.raises(exception_type) as error: + model.add_art_from_files(outcar_files, "outcar") + assert in_reason in str(error.value) diff --git a/test/tests/test_spectrum.py b/test/tests/test_spectrum.py index 125ba48..0e16126 100644 --- a/test/tests/test_spectrum.py +++ b/test/tests/test_spectrum.py @@ -9,6 +9,7 @@ import pytest from ramannoodle.polarizability.interpolation import InterpolationModel +from ramannoodle.polarizability.art import ARTModel from ramannoodle.symmetry import StructuralSymmetry from ramannoodle import io from ramannoodle.spectrum.spectrum_utils import ( @@ -64,12 +65,12 @@ def _validate_polarizabilities(model: InterpolationModel, data_directory: str) - ], indirect=["outcar_symmetry_fixture"], ) -def test_spectrum( +def test_interpolation_spectrum( outcar_symmetry_fixture: StructuralSymmetry, data_directory: str, dof_eps_outcars: list[str], ) -> None: - """Test a full spectrum calculation.""" + """Test a full spectrum calculation using InterpolationModel.""" # Setup model symmetry = outcar_symmetry_fixture _, polarizability = io.read_positions_and_polarizability( @@ -105,6 +106,63 @@ def test_spectrum( assert np.isclose(intensities, known_intensities, atol=1e-4).all() +@pytest.mark.parametrize( + "outcar_symmetry_fixture,data_directory,dof_eps_outcars", + [ + ( + "test/data/TiO2/phonons_OUTCAR", + "test/data/TiO2/", + [ + ["Ti5_0.1z_eps_OUTCAR"], + ["Ti5_0.1x_eps_OUTCAR"], + [ + "O43_0.1z_eps_OUTCAR", + "O43_m0.1z_eps_OUTCAR", + ], + ["O43_0.1x_eps_OUTCAR"], + ["O43_0.1y_eps_OUTCAR"], + ], + ), + ], + indirect=["outcar_symmetry_fixture"], +) +def test_art_spectrum( + outcar_symmetry_fixture: StructuralSymmetry, + data_directory: str, + dof_eps_outcars: list[str], +) -> None: + """Test a full spectrum calculation using ARTModel.""" + # Setup model + symmetry = outcar_symmetry_fixture + _, polarizability = io.read_positions_and_polarizability( + f"{data_directory}/ref_eps_OUTCAR", file_format="outcar" + ) + model = ARTModel(symmetry, polarizability) + for outcar_names in dof_eps_outcars: + model.add_art_from_files( + [f"{data_directory}/{name}" for name in outcar_names], file_format="outcar" + ) + + # Spectrum test + with np.load(f"{data_directory}/known_art_spectrum.npz") as known_spectrum: + phonons = io.read_phonons( + f"{data_directory}/phonons_OUTCAR", file_format="outcar" + ) + spectrum = phonons.get_raman_spectrum(model) + wavenumbers, intensities = spectrum.measure( + laser_correction=True, + laser_wavelength=532, + bose_einstein_correction=True, + temperature=300, + ) + + known_wavenumbers = known_spectrum["wavenumbers"] + known_intensities = known_spectrum["intensities"] + + assert np.isclose(wavenumbers, known_wavenumbers).all() + assert np.isclose(intensities, known_intensities, atol=1e-4).all() + + @pytest.mark.parametrize( "spectrum_path,known_gaussian_spectrum_path,known_lorentzian_spectrum_path", [