Line data    Source code

       1             : # SPDX-FileCopyrightText: 2024 PairInteraction Developers
       2             : # SPDX-License-Identifier: LGPL-3.0-or-later
       3             : 
       4           1 : from __future__ import annotations
       5             : 
       6           1 : from typing import TYPE_CHECKING
       7             : 
       8           1 : import numpy as np
       9           1 : import scipy.constants as const
      10           1 : from typing_extensions import override
      11             : 
      12           1 : from pairinteraction.green_tensor.dynamic_green_tensor import dynamic_green_tensor_homogeneous
      13           1 : from pairinteraction.green_tensor.green_tensor_base import GreenTensorBase
      14           1 : from pairinteraction.green_tensor.utils import evaluate_relative_permittivity
      15           1 : from pairinteraction.units import ureg
      16             : 
      17             : if TYPE_CHECKING:
      18             :     from typing_extensions import Self
      19             : 
      20             :     from pairinteraction.green_tensor.utils import PermittivityLike
      21             :     from pairinteraction.units import ArrayLike, NDArray, PintArrayLike
      22             : 
      23             : 
      24           1 : class GreenTensorFreeSpace(GreenTensorBase):
      25             :     """Green tensor for two atoms in free space.
      26             : 
      27             :     Examples:
      28             :         >>> from pairinteraction.green_tensor import GreenTensorFreeSpace
      29             :         >>> gt = GreenTensorFreeSpace([0, 0, 0], [10, 0, 0], unit="micrometer")
      30             :         >>> transition_energy = 2  # h * GHz
      31             :         >>> gt_dipole_dipole = gt.get(1, 1, transition_energy, "planck_constant * GHz")
      32             :         >>> print(f"{gt_dipole_dipole[0, 0]:.2f}")
      33             :         -4.79 / bohr
      34             : 
      35             :     """
      36             : 
      37           1 :     def __init__(
      38             :         self,
      39             :         pos1: ArrayLike | PintArrayLike,
      40             :         pos2: ArrayLike | PintArrayLike,
      41             :         unit: str | None = None,
      42             :         static_limit: bool = True,
      43             :         interaction_order: int = 3,
      44             :         *,
      45             :         without_vacuum_contribution: bool = False,
      46             :     ) -> None:
      47             :         """Create a Green tensor for two atoms in free space.
      48             : 
      49             :         Args:
      50             :             pos1: Position of the first atom in the given unit.
      51             :             pos2: Position of the second atom in the given unit.
      52             :             unit: The unit of the distance, e.g. "micrometer".
      53             :                 Default None expects a `pint.Quantity`.
      54             :             static_limit: If True, the static limit is used.
      55             :                 Default True.
      56             :             interaction_order: The order of interaction, e.g., 3 for dipole-dipole.
      57             :                 Defaults to 3.
      58             :             without_vacuum_contribution: Reserved for API consistency. Free-space Green tensors
      59             :                 always include the vacuum contribution, so setting this to True is invalid.
      60             : 
      61             :         """
      62           1 :         super().__init__(
      63             :             pos1, pos2, unit, static_limit, interaction_order, without_vacuum_contribution=without_vacuum_contribution
      64             :         )
      65             : 
      66           1 :     def set_relative_permittivity(self, epsilon: PermittivityLike) -> Self:
      67             :         """Set the relative permittivity of the system.
      68             : 
      69             :         Args:
      70             :             epsilon: The relative permittivity (dimensionless) of the free space. Default is 1.
      71             : 
      72             : 
      73             :         """
      74           1 :         self.epsilon = epsilon
      75           1 :         return self
      76             : 
      77           1 :     @override
      78           1 :     def _get_scaled_au(self, kappa1: int, kappa2: int, transition_energy_au: float) -> NDArray:
      79           1 :         if kappa1 == 1 and kappa2 == 1:
      80           1 :             return self._get_scaled_dipole_dipole_au(transition_energy_au)
      81           0 :         raise NotImplementedError("Only dipole-dipole Green tensors are currently implemented.")
      82             : 
      83           1 :     def _get_scaled_dipole_dipole_au(self, transition_energy_au: float) -> NDArray:
      84             :         """Calculate the dipole dipole Green tensor in cartesian coordinates for free space in atomic units.
      85             : 
      86             :         Args:
      87             :             transition_energy_au: The transition energy in atomic units at which to evaluate the Green tensor.
      88             : 
      89             :         Returns:
      90             :             The dipole dipole Green tensor in cartesian coordinates as a 3x3 array in atomic units (i.e. 1/bohr).
      91             : 
      92             :         """
      93           1 :         au_to_meter: float = ureg.Quantity(1, "atomic_unit_of_length").to("meter").magnitude
      94           1 :         pos1_m = np.array(self.pos1_au) * au_to_meter
      95           1 :         pos2_m = np.array(self.pos2_au) * au_to_meter
      96           1 :         epsilon = evaluate_relative_permittivity(self.epsilon, transition_energy_au, "hartree")
      97             : 
      98           1 :         omega_hz = ureg.Quantity(transition_energy_au, "hartree").to("hbar Hz", "spectroscopy").magnitude
      99             : 
     100             :         # unit: # m^(-3) [hbar]^(-1) [epsilon_0]^(-1)
     101           1 :         if self.without_vacuum_contribution:
     102           0 :             raise ValueError("The Green tensor for free space cannot be provided without the vacuum contribution.")
     103           1 :         gt = dynamic_green_tensor_homogeneous(pos1_m, pos2_m, omega_hz, epsilon, only_real_part=True)
     104           1 :         to_au = au_to_meter ** (-3) * ((4 * np.pi) ** (-1)) / (const.epsilon_0 * const.hbar)
     105             :         # hbar = 1, epsilon_0 = (4*np.pi)**(-1) in atomic units
     106           1 :         return np.real(gt) / to_au