From 6c19df39bb7fdb90b86f466b6ce5ae234c378e84 Mon Sep 17 00:00:00 2001 From: johannes-moegerle Date: Tue, 10 Mar 2026 11:59:39 +0100 Subject: [PATCH 01/18] add RydbergRitzParameters TypeAlias and upadte species object --- src/rydstate/species/species_object.py | 6 +++++- .../species/sqdt/species_object_sqdt.py | 20 +++++++++---------- src/rydstate/species/utils.py | 5 +++++ 3 files changed, 20 insertions(+), 11 deletions(-) diff --git a/src/rydstate/species/species_object.py b/src/rydstate/species/species_object.py index e44e05a..d4c1b42 100644 --- a/src/rydstate/species/species_object.py +++ b/src/rydstate/species/species_object.py @@ -2,7 +2,7 @@ import inspect import logging -from abc import ABC +from abc import ABC, abstractmethod from functools import cache, cached_property from typing import TYPE_CHECKING, ClassVar, TypeVar, overload @@ -125,6 +125,10 @@ def get_available_species(cls) -> list[str]: """ return sorted([subclass.name for subclass in cls._get_concrete_subclasses()]) + @property + @abstractmethod + def reference_ionization_energy_au(self) -> float: ... + @overload def get_corrected_rydberg_constant(self, unit: None = None) -> PintFloat: ... diff --git a/src/rydstate/species/sqdt/species_object_sqdt.py b/src/rydstate/species/sqdt/species_object_sqdt.py index 2e9bc9a..af15d88 100644 --- a/src/rydstate/species/sqdt/species_object_sqdt.py +++ b/src/rydstate/species/sqdt/species_object_sqdt.py @@ -17,13 +17,14 @@ from pathlib import Path from rydstate.angular.angular_ket import AngularKetBase + from rydstate.species.utils import RydbergRitzParameters from rydstate.units import PintFloat logger = logging.getLogger(__name__) class SpeciesObjectSQDT(SpeciesObject): - """Abstract base class for all species objects. + """Abstract base class for all sqdt species objects. For the electronic ground state configurations and sorted shells, see e.g. https://www.webelements.com/atoms.html @@ -41,7 +42,7 @@ class SpeciesObjectSQDT(SpeciesObject): """Ionization energy with uncertainty and unit: (value, uncertainty, unit).""" # Parameters for the extended Rydberg Ritz formula, see calc_nu - _quantum_defects: ClassVar[dict[tuple[int, float, float], tuple[float, float, float, float, float]] | None] = None + _quantum_defects: ClassVar[dict[tuple[int, float, float], RydbergRitzParameters] | None] = None """Dictionary containing the quantum defects for each (l, j_tot, s_tot) combination, i.e. _quantum_defects[(l,j_tot,s_tot)] = (d0, d2, d4, d6, d8) """ @@ -126,9 +127,6 @@ def _setup_nist_energy_levels(self, file: Path) -> None: # noqa: C901, PLR0912 if len(self._nist_energy_levels) == 0: raise ValueError(f"No NIST energy levels found for species {self.name} in file {file}.") - def __repr__(self) -> str: - return f"{self.__class__.__name__}()" - def is_allowed_shell(self, n: int, l: int, s_tot: float | None = None) -> bool: """Check if the quantum numbers describe an allowed shell. @@ -183,7 +181,7 @@ def get_ionization_energy(self, unit: str | None = None) -> PintFloat | float: return ionization_energy.to(unit, "spectroscopy").magnitude @cached_property - def ionization_energy_au(self) -> float: + def reference_ionization_energy_au(self) -> float: """Ionization energy in atomic units (Hartree).""" return self.get_ionization_energy("hartree") @@ -231,7 +229,9 @@ def calc_nu( if n <= nist_n_max and use_nist_data: # try to use NIST data if (n, l, j_tot, s_tot) in self._nist_energy_levels: energy_au = self._nist_energy_levels[(n, l, j_tot, s_tot)] - energy_au -= self.ionization_energy_au # use the cached ionization energy for better performance + energy_au -= ( + self.reference_ionization_energy_au + ) # use the cached ionization energy for better performance return calc_nu_from_energy(self.reduced_mass_au, energy_au) logger.debug( "NIST energy levels for (n=%d, l=%d, j_tot=%s, s_tot=%s) not found, using quantum defect theory.", @@ -240,10 +240,10 @@ def calc_nu( if self._quantum_defects is None: raise ValueError(f"No quantum defect data available for species {self.name}.") - quantum_defects = list(self._quantum_defects.get((l, j_tot, s_tot), [])) - if len(quantum_defects) > 5: + quantum_defects = self._quantum_defects.get((l, j_tot, s_tot), []) + if isinstance(quantum_defects, float) or np.isscalar(quantum_defects) or len(quantum_defects) > 5: raise ValueError(f"Quantum defects for {self.name} must be a list with up to 5 elements.") - d0, d2, d4, d6, d8 = quantum_defects + [0] * (5 - len(quantum_defects)) + d0, d2, d4, d6, d8 = list(quantum_defects) + [0] * (5 - len(quantum_defects)) delta_nlj = d0 + d2 / (n - d0) ** 2 + d4 / (n - d0) ** 4 + d6 / (n - d0) ** 6 + d8 / (n - d0) ** 8 return n - delta_nlj diff --git a/src/rydstate/species/utils.py b/src/rydstate/species/utils.py index f92fe98..0706b5d 100644 --- a/src/rydstate/species/utils.py +++ b/src/rydstate/species/utils.py @@ -1,5 +1,10 @@ +from __future__ import annotations + import math import re +from typing import TypeAlias + +RydbergRitzParameters: TypeAlias = tuple[float, float, float, float, float] | list[float] | float def calc_nu_from_energy(reduced_mass_au: float, energy_au: float) -> float: From df04e916d474702edbf39f21beac55e7ee5c41eb Mon Sep 17 00:00:00 2001 From: johannes-moegerle Date: Tue, 10 Mar 2026 12:09:00 +0100 Subject: [PATCH 02/18] add angular ket dummy and unkown qn --- src/rydstate/angular/angular_ket.py | 10 ++ src/rydstate/angular/angular_ket_dummy.py | 120 ++++++++++++++++++++++ src/rydstate/angular/angular_state.py | 13 ++- src/rydstate/angular/utils.py | 81 ++++++++++++--- 4 files changed, 207 insertions(+), 17 deletions(-) create mode 100644 src/rydstate/angular/angular_ket_dummy.py diff --git a/src/rydstate/angular/angular_ket.py b/src/rydstate/angular/angular_ket.py index a0ca7ad..b3c9131 100644 --- a/src/rydstate/angular/angular_ket.py +++ b/src/rydstate/angular/angular_ket.py @@ -17,7 +17,9 @@ NotSet, check_spin_addition_rule, get_possible_quantum_number_values, + is_dummy_ket, is_not_set, + is_unknown, minus_one_pow, try_trivial_spin_addition, ) @@ -232,6 +234,11 @@ def to_state(self, coupling_scheme: CouplingScheme | None = None) -> AngularStat The angular state in the specified coupling scheme. """ + if any(is_unknown(qn) for qn in self.quantum_numbers): + from rydstate.angular.angular_ket_dummy import AngularKetDummy # noqa: PLC0415 + + return AngularKetDummy(str(self), f_tot=self.f_tot, m=self.m).to_state(coupling_scheme) + if coupling_scheme is None or coupling_scheme == self.coupling_scheme: return self._create_angular_state([1], [self]) if coupling_scheme == "LS": @@ -355,6 +362,9 @@ def calc_reduced_overlap(self, other: AngularKetBase) -> float: If the kets are of different types, the overlap is calculated using the corresponding Clebsch-Gordan coefficients (/ Wigner-j symbols). """ + if any(is_dummy_ket(s) for s in [self, other]): + return 1.0 if self == other else 0.0 + for q in set(self.quantum_number_names) & set(other.quantum_number_names): if self.get_qn(q) != other.get_qn(q): return 0 diff --git a/src/rydstate/angular/angular_ket_dummy.py b/src/rydstate/angular/angular_ket_dummy.py new file mode 100644 index 0000000..07775cd --- /dev/null +++ b/src/rydstate/angular/angular_ket_dummy.py @@ -0,0 +1,120 @@ +from __future__ import annotations + +import logging +from typing import TYPE_CHECKING, ClassVar + +from rydstate.angular.angular_ket import AngularKetBase +from rydstate.angular.angular_matrix_element import is_angular_operator_type +from rydstate.angular.core_ket_base import CoreKetDummy +from rydstate.angular.utils import InvalidQuantumNumbersError, NotSet, Unknown, is_not_set + +if TYPE_CHECKING: + from typing_extensions import Self + + from rydstate.angular.angular_matrix_element import AngularMomentumQuantumNumbers, AngularOperatorType + from rydstate.angular.angular_state import AngularState + from rydstate.angular.utils import CouplingScheme, UnknownType + +logger = logging.getLogger(__name__) + + +class AngularKetDummy(AngularKetBase): + """Dummy spin ket for unknown quantum numbers.""" + + __slots__ = ("name",) + quantum_number_names: ClassVar = ("f_tot",) + coupled_quantum_numbers: ClassVar = {} + coupling_scheme = "Dummy" # type: ignore [assignment] + + name: str + """Name of the dummy ket.""" + + def __init__( + self, + name: str, + f_tot: float, + m: float | NotSet = NotSet, + ) -> None: + """Initialize the Spin ket.""" + self.name = name + + self.f_tot = f_tot + self.m = m + + super()._post_init() + + def sanity_check(self, msgs: list[str] | None = None) -> None: + """Check that the quantum numbers are valid.""" + msgs = msgs if msgs is not None else [] + + if not is_not_set(self.m) and not -self.f_tot <= self.m <= self.f_tot: + msgs.append(f"m must be between -f_tot and f_tot, but {self.f_tot=}, {self.m=}") + + if msgs: + msg = "\n ".join(msgs) + raise InvalidQuantumNumbersError(self, msg) + + def __repr__(self) -> str: + args = f"{self.name}, f_tot={self.f_tot}" + if not is_not_set(self.m): + args += f", m={self.m}" + return f"{self.__class__.__name__}({args})" + + def __str__(self) -> str: + return self.__repr__().replace("AngularKet", "") + + def __eq__(self, other: object) -> bool: + if not isinstance(other, AngularKetBase): + return NotImplemented + if not isinstance(other, AngularKetDummy): + return False + return self.name == other.name and self.f_tot == other.f_tot and self.m == other.m + + def __hash__(self) -> int: + return hash((self.name, self.f_tot, self.m)) + + def calc_reduced_matrix_element( + self: Self, + other: AngularKetBase, # noqa: ARG002 + operator: AngularOperatorType, + kappa: int, # noqa: ARG002 + ) -> float: + if not is_angular_operator_type(operator): + raise NotImplementedError(f"calc_reduced_matrix_element is not implemented for operator {operator}.") + + # ignore contributions from dummy kets + return 0 + + def get_core_ket(self) -> CoreKetDummy: + """Get the dummy core ket corresponding to this Dummy ket.""" + core_name = self.name.split("nl")[0] + return CoreKetDummy(name=core_name) + + def calc_exp_qn(self, qn: AngularMomentumQuantumNumbers) -> UnknownType: + if qn == "f_tot": + return self.f_tot + return Unknown + + def calc_std_qn(self, qn: AngularMomentumQuantumNumbers) -> UnknownType: + if qn == "f_tot": + return 0 + return Unknown + + def get_qn(self, qn: AngularMomentumQuantumNumbers) -> UnknownType: + if qn == "f_tot": + return self.f_tot + return Unknown + + def to_state(self, coupling_scheme: CouplingScheme | None = None) -> AngularState[AngularKetDummy]: # type: ignore [override] + """Convert to state in the specified coupling scheme. + + Args: + coupling_scheme: The coupling scheme to convert to (e.g. "LS", "JJ", "FJ"). + If None, the state will be a trivial state (one component) in the current coupling scheme. + + Returns: + The angular state in the specified coupling scheme. + + """ + coupling_scheme # noqa: B018 + return self._create_angular_state([1], [self]) diff --git a/src/rydstate/angular/angular_state.py b/src/rydstate/angular/angular_state.py index 7bd5ff8..cdc503c 100644 --- a/src/rydstate/angular/angular_state.py +++ b/src/rydstate/angular/angular_state.py @@ -13,13 +13,14 @@ AngularKetLS, ) from rydstate.angular.angular_matrix_element import is_angular_momentum_quantum_number -from rydstate.angular.utils import is_not_set +from rydstate.angular.utils import is_dummy_ket, is_not_set, is_unknown if TYPE_CHECKING: from collections.abc import Iterator, Sequence from typing_extensions import Self + from rydstate.angular.angular_ket_dummy import AngularKetDummy from rydstate.angular.angular_matrix_element import AngularMomentumQuantumNumbers, AngularOperatorType from rydstate.angular.utils import CouplingScheme from rydstate.units import NDArray @@ -34,7 +35,7 @@ class AngularState(Generic[_AngularKet]): def __init__( self, coefficients: Sequence[float] | NDArray, - kets: Sequence[_AngularKet], + kets: Sequence[_AngularKet | AngularKetDummy], *, warn_if_not_normalized: bool = True, ) -> None: @@ -45,7 +46,7 @@ def __init__( raise ValueError("Length of coefficients and kets must be the same.") if len(kets) == 0: raise ValueError("At least one ket must be provided.") - if not all(type(ket) is type(kets[0]) for ket in kets): + if not all((type(ket) is type(kets[0])) or is_dummy_ket(ket) for ket in kets): raise ValueError("All kets must be of the same type.") if len(set(kets)) != len(kets): raise ValueError("AngularState initialized with duplicate kets.") @@ -54,7 +55,7 @@ def __init__( if self.norm > 1: self.coefficients /= self.norm - def __iter__(self) -> Iterator[tuple[float, _AngularKet]]: + def __iter__(self) -> Iterator[tuple[float, _AngularKet | AngularKetDummy]]: return zip(self.coefficients, self.kets, strict=True).__iter__() def __repr__(self) -> str: @@ -68,7 +69,7 @@ def __str__(self) -> str: @property def coupling_scheme(self) -> CouplingScheme: """Return the coupling scheme of the state.""" - return self.kets[0].coupling_scheme + return next(ket.coupling_scheme for ket in self.kets if not is_dummy_ket(ket)) @property def norm(self) -> float: @@ -121,6 +122,7 @@ def calc_exp_qn(self, q: AngularMomentumQuantumNumbers) -> float: qs = np.array([ket.get_qn(q) for ket in self.kets]) if all(q_val == qs[0] for q_val in qs): return qs[0] # type: ignore [no-any-return] + qs = np.array([q if not is_unknown(q) else 0 for q in qs]) return np.sum(np.conjugate(self.coefficients) * self.coefficients * qs) # type: ignore [no-any-return] @@ -139,6 +141,7 @@ def calc_std_qn(self, q: AngularMomentumQuantumNumbers) -> float: qs = np.array([ket.get_qn(q) for ket in self.kets]) if all(q_val == qs[0] for q_val in qs): return 0 + qs = np.array([q if not is_unknown(q) else 0 for q in qs]) coefficients2 = np.conjugate(self.coefficients) * self.coefficients exp_q = np.sum(coefficients2 * qs) diff --git a/src/rydstate/angular/utils.py b/src/rydstate/angular/utils.py index 10382c9..0742366 100644 --- a/src/rydstate/angular/utils.py +++ b/src/rydstate/angular/utils.py @@ -1,9 +1,10 @@ from __future__ import annotations import contextlib +import logging import typing as t from functools import lru_cache -from typing import TYPE_CHECKING, Any, Literal +from typing import TYPE_CHECKING, Any, Literal, TypeAlias import numpy as np @@ -11,8 +12,12 @@ from typing_extensions import TypeIs from rydstate.angular.angular_ket import AngularKetBase + from rydstate.angular.angular_ket_dummy import AngularKetDummy from rydstate.species.species_object import SpeciesObject +logger = logging.getLogger(__name__) + + CouplingScheme = Literal["LS", "JJ", "FJ"] @@ -29,13 +34,31 @@ class NotSet(t.Protocol): def __not_set() -> None: ... +UnknownType: TypeAlias = float # type(np.nan) == float +Unknown = np.nan + + def is_not_set(obj: Any) -> TypeIs[NotSet]: # noqa: ANN401 """Check if the obj is the NotSet singleton.""" return obj is NotSet +def is_unknown(x: float | UnknownType | None) -> TypeIs[UnknownType]: + """Check if x is Unknown.""" + if x is None: + return False + return bool(np.isnan(x)) + + +def is_dummy_ket(ket: AngularKetBase) -> TypeIs[AngularKetDummy]: + """Check if ket is a AngularKetDummy.""" + from rydstate.angular.angular_ket_dummy import AngularKetDummy # noqa: PLC0415 + + return isinstance(ket, AngularKetDummy) + + class InvalidQuantumNumbersError(ValueError): - def __init__(self, ket: AngularKetBase, msg: str = "") -> None: + def __init__(self, ket: AngularKetBase | None, msg: str = "") -> None: _msg = f"Invalid quantum numbers for {ket!r}" if len(msg) > 0: _msg += f"\n {msg}" @@ -51,28 +74,55 @@ def minus_one_pow(n: float) -> int: raise ValueError(f"minus_one_pow: Invalid input {n=} is not an integer.") -def try_trivial_spin_addition(s_1: float, s_2: float, s_tot: float | None, name: str) -> float: +def try_trivial_spin_addition( + s_1: float | UnknownType, s_2: float | UnknownType, s_tot: float | UnknownType | None, name: str +) -> float | UnknownType: """Try to determine s_tot from s_1 and s_2 if it is not given. + If s_tot is Unknown, return Unknown. + If s_tot is None and any part is Unknown, return Unknown. If s_tot is None and cannot be uniquely determined from s_1 and s_2, raise an error. Otherwise return s_tot or the trivial sum s_1 + s_2. """ - if s_tot is None: - if s_1 != 0 and s_2 != 0: - msg = f"{name} must be set if both parts ({s_1} and {s_2}) are non-zero." - raise ValueError(msg) - s_tot = s_1 + s_2 - return float(s_tot) - - -def check_spin_addition_rule(s_1: float, s_2: float, s_tot: float) -> bool: + if is_unknown(s_1) or is_unknown(s_2): + if s_tot is not None: + return s_tot + return Unknown + + if s_1 != 0 and s_2 != 0: + if s_tot is not None: + return s_tot + msg = f"{name} must be set if both parts ({s_1} and {s_2}) are non-zero." + raise ValueError(msg) + + if s_1 == 0 or s_2 == 0: + calculated_s_tot = s_1 + s_2 + if is_unknown(s_tot): + logger.warning( + "%s is Unknown but could be uniquely determined as %s, using this value now", name, calculated_s_tot + ) + elif s_tot is None: + pass + elif s_tot != calculated_s_tot: + msg = f"{name} is {s_tot} but should be {calculated_s_tot} for the given parts ({s_1} and {s_2})." + raise InvalidQuantumNumbersError(None, msg) + return calculated_s_tot + + raise RuntimeError("This should never happen, all cases should be covered by the above conditions.") + + +def check_spin_addition_rule(s_1: float | UnknownType, s_2: float | UnknownType, s_tot: float | UnknownType) -> bool: r"""Check if the spin addition rule is satisfied. This means check the following conditions: :math:`|s_1 - s_2| \leq s_{tot} \leq s_1 + s_2` and :math:`s_1 + s_2 + s_{tot}` is an integer + + If any of the quantum numbers is Unknown, return True (cannot verify). """ + if is_unknown(s_1) or is_unknown(s_2) or is_unknown(s_tot): + return True return abs(s_1 - s_2) <= s_tot <= s_1 + s_2 and (s_1 + s_2 + s_tot) % 1 == 0 @@ -83,6 +133,13 @@ def get_possible_quantum_number_values(s_1: float, s_2: float, s_tot: float | No return [float(s) for s in np.arange(abs(s_1 - s_2), s_1 + s_2 + 1, 1)] +def is_equal(qn1: float | UnknownType, qn2: float | UnknownType) -> bool: + """Check if two quantum numbers are equal, treating Unknown equal to Unknown.""" + if is_unknown(qn1) and is_unknown(qn2): + return True + return qn1 == qn2 + + @lru_cache(maxsize=1_000) def quantum_numbers_to_angular_ket( species: str | SpeciesObject, From 1a721bb19e773c4a89574c342c852f0175ef0cca Mon Sep 17 00:00:00 2001 From: johannes-moegerle Date: Tue, 10 Mar 2026 12:09:38 +0100 Subject: [PATCH 03/18] add angular core ket --- src/rydstate/angular/angular_ket.py | 5 ++ src/rydstate/angular/core_ket_base.py | 71 +++++++++++++++++++++++++++ 2 files changed, 76 insertions(+) create mode 100644 src/rydstate/angular/core_ket_base.py diff --git a/src/rydstate/angular/angular_ket.py b/src/rydstate/angular/angular_ket.py index b3c9131..74ecb58 100644 --- a/src/rydstate/angular/angular_ket.py +++ b/src/rydstate/angular/angular_ket.py @@ -12,6 +12,7 @@ is_angular_momentum_quantum_number, is_angular_operator_type, ) +from rydstate.angular.core_ket_base import CoreKet from rydstate.angular.utils import ( InvalidQuantumNumbersError, NotSet, @@ -740,3 +741,7 @@ def sanity_check(self, msgs: list[str] | None = None) -> None: msgs.append(f"{self.f_c=}, {self.j_r=}, {self.f_tot=} don't satisfy spin addition rule.") super().sanity_check(msgs) + + def get_core_ket(self) -> CoreKet: + """Get the core ket corresponding to this FJ ket.""" + return CoreKet(i_c=self.i_c, s_c=self.s_c, l_c=self.l_c, j_c=self.j_c, f_c=self.f_c) diff --git a/src/rydstate/angular/core_ket_base.py b/src/rydstate/angular/core_ket_base.py new file mode 100644 index 0000000..0d3f0a2 --- /dev/null +++ b/src/rydstate/angular/core_ket_base.py @@ -0,0 +1,71 @@ +from __future__ import annotations + +from typing import TYPE_CHECKING + +from rydstate.angular.utils import is_equal, is_unknown, try_trivial_spin_addition + +if TYPE_CHECKING: + from rydstate.angular.utils import UnknownType + + +class CoreKet: + __slots__ = ("i_c", "s_c", "l_c", "j_c", "f_c") + + def __init__( + self, + i_c: float | None = None, + s_c: float | None = None, + l_c: int | UnknownType | None = None, + j_c: float | None = None, + f_c: float | None = None, + ) -> None: + """Initialize the core angular ket.""" + if i_c is None: + raise ValueError("Nuclear spin i_c must be set.") + self.i_c = float(i_c) + + if s_c is None: + raise ValueError("Core spin s_c must be set.") + self.s_c = float(s_c) + + if l_c is None: + raise ValueError("Core orbital angular momentum l_c must be set.") + self.l_c = int(l_c) if not is_unknown(l_c) else l_c + + self.j_c = try_trivial_spin_addition(self.l_c, self.s_c, j_c, "j_c") + self.f_c = try_trivial_spin_addition(self.j_c, self.i_c, f_c, "f_c") + + def __repr__(self) -> str: + return f"CoreKet(i_c={self.i_c}, s_c={self.s_c}, l_c={self.l_c}, j_c={self.j_c}, f_c={self.f_c})" + + def __hash__(self) -> int: + return hash((str(type(self)), self.i_c, self.s_c, self.l_c, self.j_c, self.f_c)) + + def __eq__(self, other: object) -> bool: + if not isinstance(other, CoreKet): + raise NotImplementedError(f"Cannot compare {self!r} with {other!r}.") + return ( + is_equal(self.i_c, other.i_c) + and is_equal(self.s_c, other.s_c) + and is_equal(self.l_c, other.l_c) + and is_equal(self.j_c, other.j_c) + and is_equal(self.f_c, other.f_c) + ) + + +class CoreKetDummy(CoreKet): + """Dummy core spin ket for unknown quantum numbers.""" + + def __init__(self, name: str) -> None: + self.name = name + + def __repr__(self) -> str: + return f"CoreKetDummy(name={self.name})" + + def __hash__(self) -> int: + return hash((str(type(self)), self.name)) + + def __eq__(self, other: object) -> bool: + if not isinstance(other, CoreKetDummy): + raise NotImplementedError(f"Cannot compare {self!r} with {other!r}.") + return self.name == other.name From b7ebbcda908de7016ce14d61dc2200ba63f83572 Mon Sep 17 00:00:00 2001 From: johannes-moegerle Date: Tue, 10 Mar 2026 12:00:09 +0100 Subject: [PATCH 04/18] add mqdt species object and parameters --- pyproject.toml | 1 + src/rydstate/species/__init__.py | 18 + src/rydstate/species/mqdt/__init__.py | 18 + src/rydstate/species/mqdt/fmodel.py | 313 +++++ .../species/mqdt/species_object_mqdt.py | 82 ++ src/rydstate/species/mqdt/sr87.py | 488 ++++++++ src/rydstate/species/mqdt/sr88.py | 303 +++++ src/rydstate/species/mqdt/yb171.py | 1043 +++++++++++++++++ src/rydstate/species/mqdt/yb173.py | 394 +++++++ src/rydstate/species/mqdt/yb174.py | 652 +++++++++++ 10 files changed, 3312 insertions(+) create mode 100644 src/rydstate/species/mqdt/__init__.py create mode 100644 src/rydstate/species/mqdt/fmodel.py create mode 100644 src/rydstate/species/mqdt/species_object_mqdt.py create mode 100644 src/rydstate/species/mqdt/sr87.py create mode 100644 src/rydstate/species/mqdt/sr88.py create mode 100644 src/rydstate/species/mqdt/yb171.py create mode 100644 src/rydstate/species/mqdt/yb173.py create mode 100644 src/rydstate/species/mqdt/yb174.py diff --git a/pyproject.toml b/pyproject.toml index 0bf8001..3d17659 100644 --- a/pyproject.toml +++ b/pyproject.toml @@ -137,6 +137,7 @@ mccabe.max-complexity = 10 [tool.ruff.lint.per-file-ignores] "docs/examples/*.ipynb" = ["T201"] +"src/rydstate/species/mqdt/*.py" = ["RUF012", "N801"] [tool.ruff.lint.isort] combine-as-imports = true diff --git a/src/rydstate/species/__init__.py b/src/rydstate/species/__init__.py index 6a96fb4..432abb7 100644 --- a/src/rydstate/species/__init__.py +++ b/src/rydstate/species/__init__.py @@ -1,3 +1,13 @@ +from rydstate.species.mqdt import ( + FModel, + FModelSQDT, + SpeciesObjectMQDT, + Strontium87MQDT, + Strontium88MQDT, + Ytterbium171MQDT, + Ytterbium173MQDT, + Ytterbium174MQDT, +) from rydstate.species.species_object import SpeciesObject from rydstate.species.sqdt import ( Cesium, @@ -17,6 +27,8 @@ __all__ = [ "Cesium", + "FModel", + "FModelSQDT", "Hydrogen", "HydrogenTextBook", "Lithium", @@ -24,10 +36,16 @@ "Rubidium", "Sodium", "SpeciesObject", + "SpeciesObjectMQDT", "SpeciesObjectSQDT", "Strontium87", + "Strontium87MQDT", "Strontium88", + "Strontium88MQDT", "Ytterbium171", + "Ytterbium171MQDT", "Ytterbium173", + "Ytterbium173MQDT", "Ytterbium174", + "Ytterbium174MQDT", ] diff --git a/src/rydstate/species/mqdt/__init__.py b/src/rydstate/species/mqdt/__init__.py new file mode 100644 index 0000000..e4ddb9a --- /dev/null +++ b/src/rydstate/species/mqdt/__init__.py @@ -0,0 +1,18 @@ +from rydstate.species.mqdt.fmodel import FModel, FModelSQDT +from rydstate.species.mqdt.species_object_mqdt import SpeciesObjectMQDT +from rydstate.species.mqdt.sr87 import Strontium87MQDT +from rydstate.species.mqdt.sr88 import Strontium88MQDT +from rydstate.species.mqdt.yb171 import Ytterbium171MQDT +from rydstate.species.mqdt.yb173 import Ytterbium173MQDT +from rydstate.species.mqdt.yb174 import Ytterbium174MQDT + +__all__ = [ + "FModel", + "FModelSQDT", + "SpeciesObjectMQDT", + "Strontium87MQDT", + "Strontium88MQDT", + "Ytterbium171MQDT", + "Ytterbium173MQDT", + "Ytterbium174MQDT", +] diff --git a/src/rydstate/species/mqdt/fmodel.py b/src/rydstate/species/mqdt/fmodel.py new file mode 100644 index 0000000..7dc250e --- /dev/null +++ b/src/rydstate/species/mqdt/fmodel.py @@ -0,0 +1,313 @@ +from __future__ import annotations + +import math +from functools import cached_property +from typing import TYPE_CHECKING, ClassVar, overload + +import numpy as np + +from rydstate.angular.utils import is_unknown +from rydstate.species.utils import calc_energy_from_nu, calc_nu_from_energy +from rydstate.utils.linalg import find_roots + +if TYPE_CHECKING: + from rydstate.angular.angular_ket import AngularKetBase, AngularKetFJ + from rydstate.angular.angular_ket_dummy import AngularKetDummy + from rydstate.species.mqdt.species_object_mqdt import SpeciesObjectMQDT + from rydstate.species.utils import RydbergRitzParameters + from rydstate.units import NDArray, PintFloat + + +class FModel: + """Class to store the parameters of a MQDT model for a given species.""" + + reference: ClassVar[str | None] = None + """Reference for the MQDT model, e.g., a publication doi where the model is described.""" + + species_name: ClassVar[str] + """The species for which the MQDT model is defined.""" + + name: ClassVar[str | None] + """The name of the atomic species.""" + + f_tot: ClassVar[float] + """Total angular momentum f_tot of the Rydberg state.""" + + nu_range: ClassVar[tuple[float, float]] + """Range of effective principal quantum numbers nu for which the MQDT model is valid.""" + + inner_channels: ClassVar[list[AngularKetBase]] + """List of inner channels in the MQDT model.""" + outer_channels: ClassVar[list[AngularKetFJ | AngularKetDummy]] + """List of outer channels in the MQDT model.""" + + eigen_quantum_defects: ClassVar[list[RydbergRitzParameters]] + mixing_angles: ClassVar[list[tuple[int, int, RydbergRitzParameters]]] + + manual_frame_transformation_outer_inner: ClassVar[NDArray | None] = None + """Optional manually specified frame transformation matrix Q mapping inner to outer channels. + This is mainly needed for models with unknown quantum numbers, + where the frame transformation cannot (yet) be computed from Wigner coefficients. + """ + + @property + def nu_min(self) -> float: + """Minimum nu for which the model is valid.""" + return self.nu_range[0] + + @property + def nu_max(self) -> float: + """Maximum nu for which the model is valid.""" + return self.nu_range[1] + + @cached_property + def species(self) -> SpeciesObjectMQDT: + """Return the SpeciesObjectMQDT associated with this model.""" + from rydstate.species.mqdt.species_object_mqdt import SpeciesObjectMQDT # noqa: PLC0415 + + return SpeciesObjectMQDT.from_name(self.species_name) + + @overload + def get_ionization_thresholds(self, unit: None = None) -> list[PintFloat]: ... + + @overload + def get_ionization_thresholds(self, unit: str) -> list[float]: ... + + def get_ionization_thresholds(self, unit: str | None = "hartree") -> list[PintFloat] | list[float]: + """Return the ionization thresholds for all channels. + + Args: + unit: Desired unit for the ionization thresholds. Default is atomic units "hartree". + + Returns: + List of ionization thresholds in the desired unit. + + """ + return [self.species.get_ionization_threshold(ket.get_core_ket(), unit=unit) for ket in self.outer_channels] # type: ignore [return-value] + + @cached_property # don't remove this caching without benchmarking it!!! + def ionization_thresholds_au(self) -> list[float]: + """Return the ionization thresholds for all channels in atomic units.""" + return self.get_ionization_thresholds(unit="hartree") + + def calc_channel_nuis(self, nu_ref: float) -> NDArray: + r"""Return the channel-dependent effective principal quantum numbers nui. + + The channel dependent effective principal quantum numbers nui are defined via + + .. math:: + E = I_i - \frac{Ry}{2 \nu_i^2} = I_{\text{ref}} - \frac{Ry}{nu_ref^2} + + Args: + nu_ref: Reference effective principal quantum number. + + Returns: + List of channel nui values. + + """ + thresholds = self.ionization_thresholds_au + ioniz_ref = self.species.reference_ionization_energy_au + reduced_mass_au = self.species.reduced_mass_au + energy_from_nu_ref = calc_energy_from_nu(reduced_mass_au, nu_ref) + nuis = [ + calc_nu_from_energy(reduced_mass_au, ioniz_ref + energy_from_nu_ref - threshold) for threshold in thresholds + ] + return np.array(nuis) + + def calc_k_matrix_closecoupling(self, nu_ref: float) -> NDArray: + r"""Return diagonal K-matrix in the close-coupling frame. + + Diagonal entries are tan(\pi * \mu_\alpha) where \mu_\alpha are the eigen quantum defects + evaluated at the channel-dependent effective principal quantum numbers nui. + + Args: + nu_ref: Reference effective principal quantum number. + + Returns: + Diagonal K-matrix in the close-coupling frame. + + """ + nuis = self.calc_channel_nuis(nu_ref) + defects = [ + calc_energy_dependent_quantity(nu, params) + for nu, params in zip(nuis, self.eigen_quantum_defects, strict=True) + ] + return np.diag(np.tan(np.pi * np.array(defects))) + + def calc_frame_transformation_outer_inner(self) -> NDArray: + """Return the frame transformation matrix Q mapping inner to outer channels. + + Computed from the overlaps (Wigner coefficients) between inner_channels and outer_channels. + + Returns: + Unitary transformation matrix Q (n_outer, n_inner). + + """ + if self.manual_frame_transformation_outer_inner is not None: + return self.manual_frame_transformation_outer_inner + + n = len(self.inner_channels) + u = np.zeros((n, n)) + for i, outer in enumerate(self.outer_channels): + for j, inner in enumerate(self.inner_channels): + if any(is_unknown(qn) for qn in outer.quantum_numbers) or any( + is_unknown(qn) for qn in inner.quantum_numbers + ): + raise ValueError( + "Cannot compute frame transformation for channels with unknown quantum numbers." + f"Please provide a manual frame transformation matrix for {self.name}." + ) + u[i, j] = outer.calc_reduced_overlap(inner) + return u + + @cached_property # don't remove this caching without benchmarking it!!! + def frame_transformation_outer_inner(self) -> NDArray: + """Cached version of calc_frame_transformation_outer_inner.""" + return self.calc_frame_transformation_outer_inner() + + def calc_frame_transformation_inner_closecoupling(self, nu_ref: float) -> NDArray: + """Return the frame transformation matrix R mapping close-coupling to inner channels. + + Computed as rotation matrix from the mixing angles. + Applies successive 2x2 rotations between the channels specified by mixing_angles. + + Args: + nu_ref: Reference effective principal quantum number. + + Returns: + Unitary transformation matrix R (n_inner, n_closecoupling). + + """ + n = len(self.inner_channels) + if len(self.mixing_angles) == 0: + return np.eye(n) + # Find reference channel nu for energy-dependent angles + # convention: first involved channel of first energy-dependent mixing entry + ref_nu: float | None = None + for i_idx, _j_idx, params in self.mixing_angles: + if isinstance(params, list) and len(params) > 1: + nuis = self.calc_channel_nuis(nu_ref) + ref_nu = float(nuis[i_idx]) + break + if ref_nu is None: + ref_nu = 0.0 # unused; angles are constant + rot = np.eye(n) + for i_idx, j_idx, params in self.mixing_angles: + angle = calc_energy_dependent_quantity(ref_nu, params) + r = np.eye(n) + r[i_idx, i_idx] = np.cos(angle) + r[i_idx, j_idx] = -np.sin(angle) + r[j_idx, i_idx] = np.sin(angle) + r[j_idx, j_idx] = np.cos(angle) + rot = rot @ r + return rot + + def calc_frame_transformation(self, nu_ref: float) -> NDArray: + """Return the full frame transformation U from close-coupling to outer channel frame. + + Combines the unitary frame transformation Q with the rotation matrix R. + + Args: + nu_ref: Reference effective principal quantum number. + + Returns: + Frame transformation matrix U = Q R (n_outer, n_closecoupling). + + """ + return self.frame_transformation_outer_inner @ self.calc_frame_transformation_inner_closecoupling(nu_ref) + + def calc_k_matrix(self, nu_ref: float) -> NDArray: + r"""Return the K-matrix in the collision (outer) channel frame. + + The K-matrix is defined as + + .. math:: + K = tan(\pi \mu) = U tan(\pi \mu_{\alpha}) U^T + + where U is the frame transformation matrix and \mu_{\alpha} are the eigen quantum defects. + + Args: + nu_ref: Reference effective principal quantum number. + + Returns: + K-matrix in the collision (outer) channel frame, K = tan(\pi \mu). + + """ + transform = self.calc_frame_transformation(nu_ref) + kbar = self.calc_k_matrix_closecoupling(nu_ref) + return transform @ kbar @ transform.T + + def calc_m_matrix(self, nu_ref: float) -> NDArray: + r"""Return the M-matrix in the collision (outer) channel frame. + + The M-matrix is defined as + + .. math:: + M = tan(β) + K = tan(\pi \nu) + tan(\pi \mu) + + Args: + nu_ref: Reference effective principal quantum number. + + Returns: + M-matrix in the collision (outer) channel frame, M = tan(β) + K. + + """ + kmat = self.calc_k_matrix(nu_ref) + nuis = self.calc_channel_nuis(nu_ref) + return np.array(np.diag(np.tan(np.pi * nuis)) + kmat) + + def _calc_scaled_m_matrix(self, nu_ref: float) -> NDArray: + # we scale the M-matrix by cos(pi * nui) to improve numerical stability when finding roots of det(M) = 0 + # this is especially important for states with nu close to half integer + kmat = self.calc_k_matrix(nu_ref) + nuis = self.calc_channel_nuis(nu_ref) + return np.array(np.diag(np.sin(np.pi * nuis)) + np.cos(np.pi * nuis) * kmat) + + def calc_detm_roots(self, nu_min: float, nu_max: float) -> list[float]: + """Find zeros of det(mmat(nu_ref)) in [nu_min, nu_max] using scipy root finding. + + Args: + nu_min: Lower bound of the search range. + nu_max: Upper bound of the search range. + + Returns: + List of nu_ref values where det(mmat) ≈ 0. + + """ + + def det_mmat(nu: float) -> float: + return float(np.linalg.det(self._calc_scaled_m_matrix(nu))) + + return find_roots(det_mmat, nu_min, nu_max) + + +class FModelSQDT(FModel): + def __init__(self, species_name: str, channel: AngularKetFJ) -> None: + self.species_name = species_name # type: ignore [misc] + self.name = f"SQDT l_r={channel.l_r}, j_r={channel.j_r}, f_tot={channel.f_tot}, nu > {channel.l_r + 1}" # type: ignore [misc] + self.f_tot = channel.f_tot # type: ignore [misc] + self.nu_range = (channel.l_r + 1, math.inf) # type: ignore [misc] + self.inner_channels = [channel] # type: ignore [misc] + self.outer_channels = [channel] # type: ignore [misc] + self.eigen_quantum_defects = [0] # type: ignore [misc] + self.mixing_angles = [] # type: ignore [misc] + + +def calc_energy_dependent_quantity(nui: float, qantity: RydbergRitzParameters) -> float: + """Evaluate a quantity at channel nu using polynomial convention: q₀ + q₁/ν² + q₂/ν⁴ + ... + + Args: + nui: Channel-dependent effective principal quantum number. + qantity: Quantity parameters. A single float is a constant value; a list gives + polynomial coefficients [q₀, q₁, q₂, ...]. + + Returns: + Quantity value at the given nu. + + """ + if isinstance(qantity, (int, float)): + return float(qantity) + result = 0.0 + for i, coeff in enumerate(qantity): + result += float(coeff) * (1.0 / nui**2) ** i + return result diff --git a/src/rydstate/species/mqdt/species_object_mqdt.py b/src/rydstate/species/mqdt/species_object_mqdt.py new file mode 100644 index 0000000..8592e18 --- /dev/null +++ b/src/rydstate/species/mqdt/species_object_mqdt.py @@ -0,0 +1,82 @@ +from __future__ import annotations + +from functools import cached_property +from typing import TYPE_CHECKING, ClassVar, overload + +from rydstate.angular.core_ket_base import CoreKet, CoreKetDummy +from rydstate.angular.utils import Unknown +from rydstate.species.mqdt.fmodel import FModel, FModelSQDT +from rydstate.species.species_object import SpeciesObject +from rydstate.units import ureg + +if TYPE_CHECKING: + from rydstate.angular.angular_ket import AngularKetFJ + from rydstate.units import PintFloat + + +class SpeciesObjectMQDT(SpeciesObject): + i_c: ClassVar[float] + + _ionization_threshold_dict: dict[CoreKet, tuple[float, float | None, str]] + """Dictionary containing the ionization thresholds for the different core states. + The thresholds are given in the form of a tuple (ionization_threshold, uncertainty, unit). + """ + + core_ground_state: CoreKet + """The ground state configuration of the atomic core.""" + + nuclear_dipole: float + """Nuclear dipole moment of the species.""" + + @cached_property + def models(self) -> list[FModel]: + """List of MQDT models available for the species.""" + models = [model for model in FModel.__subclasses__() if getattr(model, "species_name", None) == self.name] + return [model() for model in models] + + @overload + def get_ionization_threshold(self, core_ket: CoreKet, unit: None = None) -> PintFloat: ... + + @overload + def get_ionization_threshold(self, core_ket: CoreKet, unit: str) -> float: ... + + def get_ionization_threshold(self, core_ket: CoreKet, unit: str | None = "hartree") -> PintFloat | float: + """Return the ionization energy of the channel given by the core_ket in the desired unit. + + Args: + core_ket: The core ket for which to return the ionization energy. + unit: Desired unit for the ionization energy. Default is atomic units "hartree". + + Returns: + Ionization energy in the desired unit. + + """ + if core_ket not in self._ionization_threshold_dict: + if isinstance(core_ket, CoreKetDummy): + raise ValueError(f"Core ket {core_ket} is a dummy state. Ionization energy is not defined.") + core_ket = CoreKet(i_c=self.i_c, s_c=core_ket.s_c, l_c=core_ket.l_c, j_c=core_ket.j_c, f_c=Unknown) + if core_ket not in self._ionization_threshold_dict: + core_ket = CoreKet(i_c=self.i_c, s_c=core_ket.s_c, l_c=core_ket.l_c, j_c=Unknown, f_c=Unknown) + if core_ket not in self._ionization_threshold_dict: + raise ValueError(f"Ionization energy for core ket {core_ket} is not defined.") + + ionization_energy_tuple = self._ionization_threshold_dict[core_ket] + ionization_energy: PintFloat = ureg.Quantity(ionization_energy_tuple[0], ionization_energy_tuple[2]) + ionization_energy = ionization_energy.to("hartree", "spectroscopy") + if unit is None: + return ionization_energy + if unit == "a.u.": + return ionization_energy.magnitude + return ionization_energy.to(unit, "spectroscopy").magnitude + + @cached_property + def reference_ionization_energy_au(self) -> float: + """Ionization energy in atomic units (Hartree).""" + return self.get_ionization_threshold(self.core_ground_state, unit="hartree") + + def get_mqdt_models(self, outer_channel: AngularKetFJ) -> list[FModel]: + """Return a list of MQDT models for the outer_channel.""" + models = [model for model in self.models if any(ket == outer_channel for ket in model.outer_channels)] + if len(models) == 0: + models = [FModelSQDT(self.name, outer_channel)] + return models diff --git a/src/rydstate/species/mqdt/sr87.py b/src/rydstate/species/mqdt/sr87.py new file mode 100644 index 0000000..ea74a95 --- /dev/null +++ b/src/rydstate/species/mqdt/sr87.py @@ -0,0 +1,488 @@ +from __future__ import annotations + +import numpy as np + +from rydstate.angular.angular_ket import AngularKetFJ, AngularKetLS +from rydstate.angular.core_ket_base import CoreKet +from rydstate.species.mqdt.fmodel import FModel +from rydstate.species.mqdt.species_object_mqdt import SpeciesObjectMQDT +from rydstate.units import electron_mass, rydberg_constant + + +class Strontium87MQDT(SpeciesObjectMQDT): + name = "Sr87_mqdt" + Z = 38 + i_c = 4.5 + number_valence_electrons = 2 + + potential_type_default = "model_potential_fei_2009" + + _isotope_mass = 86.9088774970 # u + _corrected_rydberg_constant = ( + rydberg_constant.m / (1 + electron_mass.to("u").m / _isotope_mass), + None, + str(rydberg_constant.u), + ) + + _ionization_threshold_dict = { + CoreKet(i_c, 0.5, 0, 0.5, 4): (45932.287373577, None, "1/cm"), + CoreKet(i_c, 0.5, 0, 0.5, 5): (45932.120512528, None, "1/cm"), + } + core_ground_state = CoreKet(i_c, 0.5, 0, 0.5, 4) + nuclear_dipole = -1.0936030 + + +# -------------------------------------------------------- +# MQDT models valid at large n +# -------------------------------------------------------- + + +class Sr87_S35_HighN(FModel): + species_name = "Sr87_mqdt" + name = "S F=7/2, nu > 11" + f_tot = 3.5 + nu_range = (11.0, np.inf) + reference = "10.1088/1361-6455/ab3c26" + + inner_channels = [ + AngularKetLS(l_c=0, l_r=0, l_tot=0, s_tot=1, j_tot=1, f_tot=3.5, species="Sr87"), # "5sns 3S1" + ] + outer_channels = [ + AngularKetFJ(l_c=0, l_r=0, j_c=0.5, f_c=4, j_r=0.5, f_tot=3.5, species="Sr87"), + ] + + eigen_quantum_defects = [ + [3.370778, 0.418, -0.3], + ] + mixing_angles = [] + + +class Sr87_S45_HighN(FModel): + species_name = "Sr87_mqdt" + name = "S F=9/2, nu > 11" + f_tot = 4.5 + nu_range = (11.0, np.inf) + reference = "10.1088/1361-6455/ab3c26" + + inner_channels = [ + AngularKetLS(l_c=0, l_r=0, l_tot=0, s_tot=0, j_tot=0, f_tot=4.5, species="Sr87"), # "5sns 1S0" + AngularKetLS(l_c=0, l_r=0, l_tot=0, s_tot=1, j_tot=1, f_tot=4.5, species="Sr87"), # "5sns 3S1" + ] + outer_channels = [ + AngularKetFJ(l_c=0, l_r=0, j_c=0.5, f_c=4, j_r=0.5, f_tot=4.5, species="Sr87"), + AngularKetFJ(l_c=0, l_r=0, j_c=0.5, f_c=5, j_r=0.5, f_tot=4.5, species="Sr87"), + ] + + eigen_quantum_defects = [ + [3.26896, -0.138, 0.9], + [3.370778, 0.418, -0.3], + ] + mixing_angles = [] + + +class Sr87_S55_HighN(FModel): + species_name = "Sr87_mqdt" + name = "S F=11/2, nu > 11" + f_tot = 5.5 + nu_range = (11.0, np.inf) + reference = "10.1088/1361-6455/ab3c26" + + inner_channels = [ + AngularKetLS(l_c=0, l_r=0, l_tot=0, s_tot=1, j_tot=1, f_tot=5.5, species="Sr87"), # "5sns 3S1" + ] + outer_channels = [ + AngularKetFJ(l_c=0, l_r=0, j_c=0.5, f_c=5, j_r=0.5, f_tot=5.5, species="Sr87"), + ] + + eigen_quantum_defects = [ + [3.370778, 0.418, -0.3], + ] + mixing_angles = [] + + +# -------------------------------------------------------- +# Low-n models +# -------------------------------------------------------- + + +class Sr87_P45_LowN(FModel): + species_name = "Sr87_mqdt" + name = "P F=9/2 (clock), 1.8 < nu < 2.2" + f_tot = 4.5 + nu_range = (1.8, 2.2) + reference = "10.1088/1361-6455/ab3c26" + + inner_channels = [ + AngularKetLS(l_c=0, l_r=1, l_tot=1, s_tot=0, j_tot=1, f_tot=4.5, species="Sr87"), # "5snp 1P1" + AngularKetLS(l_c=0, l_r=1, l_tot=1, s_tot=1, j_tot=0, f_tot=4.5, species="Sr87"), # "5snp 3P0" + AngularKetLS(l_c=0, l_r=1, l_tot=1, s_tot=1, j_tot=1, f_tot=4.5, species="Sr87"), # "5snp 3P1" + AngularKetLS(l_c=0, l_r=1, l_tot=1, s_tot=1, j_tot=2, f_tot=4.5, species="Sr87"), # "5snp 3P2" + ] + outer_channels = [ + AngularKetFJ(l_c=0, l_r=1, j_c=0.5, f_c=4, j_r=0.5, f_tot=4.5, species="Sr87"), + AngularKetFJ(l_c=0, l_r=1, j_c=0.5, f_c=4, j_r=1.5, f_tot=4.5, species="Sr87"), + AngularKetFJ(l_c=0, l_r=1, j_c=0.5, f_c=5, j_r=0.5, f_tot=4.5, species="Sr87"), + AngularKetFJ(l_c=0, l_r=1, j_c=0.5, f_c=5, j_r=1.5, f_tot=4.5, species="Sr87"), + ] + + eigen_quantum_defects = [ + [0.8720737, 0], + [0.13689075, 0], + [0.13143188, 0], + [0.11955235, 0], + ] + mixing_angles = [] + + +# -------------------------------------------------------- +# High-n P models +# -------------------------------------------------------- + + +class Sr87_P25_HighN(FModel): + species_name = "Sr87_mqdt" + name = "P F=5/2, nu > 5" + f_tot = 2.5 + nu_range = (5.0, np.inf) + reference = "10.1088/1361-6455/ab3c26" + + inner_channels = [ + AngularKetLS(l_c=0, l_r=1, l_tot=1, s_tot=1, j_tot=2, f_tot=2.5, species="Sr87"), # "5snp 3P2" + ] + outer_channels = [ + AngularKetFJ(l_c=0, l_r=1, j_c=0.5, f_c=4, j_r=1.5, f_tot=2.5, species="Sr87"), + ] + + eigen_quantum_defects = [ + [2.882, 0.446, -1.9], + ] + mixing_angles = [] + + +class Sr87_P35_HighN(FModel): + species_name = "Sr87_mqdt" + name = "P F=7/2, nu > 5" + f_tot = 3.5 + nu_range = (5.0, np.inf) + reference = "10.1088/1361-6455/ab3c26" + + inner_channels = [ + AngularKetLS(l_c=0, l_r=1, l_tot=1, s_tot=0, j_tot=1, f_tot=3.5, species="Sr87"), # "5snp 1P1" + AngularKetLS(l_c=0, l_r=1, l_tot=1, s_tot=1, j_tot=1, f_tot=3.5, species="Sr87"), # "5snp 3P1" + AngularKetLS(l_c=0, l_r=1, l_tot=1, s_tot=1, j_tot=2, f_tot=3.5, species="Sr87"), # "5snp 3P2" + ] + outer_channels = [ + AngularKetFJ(l_c=0, l_r=1, j_c=0.5, f_c=4, j_r=0.5, f_tot=3.5, species="Sr87"), + AngularKetFJ(l_c=0, l_r=1, j_c=0.5, f_c=4, j_r=1.5, f_tot=3.5, species="Sr87"), + AngularKetFJ(l_c=0, l_r=1, j_c=0.5, f_c=5, j_r=1.5, f_tot=3.5, species="Sr87"), + ] + + eigen_quantum_defects = [ + [2.724, -4.67, -157], + [2.8826, 0.407, -1.3], + [2.882, 0.446, -1.9], + ] + mixing_angles = [] + + +class Sr87_P45_HighN(FModel): + species_name = "Sr87_mqdt" + name = "P F=9/2, nu > 7" + f_tot = 4.5 + nu_range = (7.0, np.inf) + reference = "10.1088/1361-6455/ab3c26" + + inner_channels = [ + AngularKetLS(l_c=0, l_r=1, l_tot=1, s_tot=0, j_tot=1, f_tot=4.5, species="Sr87"), # "5snp 1P1" + AngularKetLS(l_c=0, l_r=1, l_tot=1, s_tot=1, j_tot=0, f_tot=4.5, species="Sr87"), # "5snp 3P0" + AngularKetLS(l_c=0, l_r=1, l_tot=1, s_tot=1, j_tot=1, f_tot=4.5, species="Sr87"), # "5snp 3P1" + AngularKetLS(l_c=0, l_r=1, l_tot=1, s_tot=1, j_tot=2, f_tot=4.5, species="Sr87"), # "5snp 3P2" + ] + outer_channels = [ + AngularKetFJ(l_c=0, l_r=1, j_c=0.5, f_c=4, j_r=0.5, f_tot=4.5, species="Sr87"), + AngularKetFJ(l_c=0, l_r=1, j_c=0.5, f_c=4, j_r=1.5, f_tot=4.5, species="Sr87"), + AngularKetFJ(l_c=0, l_r=1, j_c=0.5, f_c=5, j_r=0.5, f_tot=4.5, species="Sr87"), + AngularKetFJ(l_c=0, l_r=1, j_c=0.5, f_c=5, j_r=1.5, f_tot=4.5, species="Sr87"), + ] + + eigen_quantum_defects = [ + [2.724, -4.67, -157], + [2.8867, 0.44, -1.9], + [2.8826, 0.407, -1.3], + [2.882, 0.446, -1.9], + ] + mixing_angles = [] + + +class Sr87_P55_HighN(FModel): + species_name = "Sr87_mqdt" + name = "P F=11/2, nu > 5" + f_tot = 5.5 + nu_range = (5.0, np.inf) + reference = "10.1088/1361-6455/ab3c26" + + inner_channels = [ + AngularKetLS(l_c=0, l_r=1, l_tot=1, s_tot=0, j_tot=1, f_tot=5.5, species="Sr87"), # "5snp 1P1" + AngularKetLS(l_c=0, l_r=1, l_tot=1, s_tot=1, j_tot=1, f_tot=5.5, species="Sr87"), # "5snp 3P1" + AngularKetLS(l_c=0, l_r=1, l_tot=1, s_tot=1, j_tot=2, f_tot=5.5, species="Sr87"), # "5snp 3P2" + ] + outer_channels = [ + AngularKetFJ(l_c=0, l_r=1, j_c=0.5, f_c=4, j_r=1.5, f_tot=5.5, species="Sr87"), + AngularKetFJ(l_c=0, l_r=1, j_c=0.5, f_c=5, j_r=0.5, f_tot=5.5, species="Sr87"), + AngularKetFJ(l_c=0, l_r=1, j_c=0.5, f_c=5, j_r=1.5, f_tot=5.5, species="Sr87"), + ] + + eigen_quantum_defects = [ + [2.724, -4.67, -157], + [2.8826, 0.407, -1.3], + [2.882, 0.446, -1.9], + ] + mixing_angles = [] + + +class Sr87_P65_HighN(FModel): + species_name = "Sr87_mqdt" + name = "P F=13/2, nu > 5" + f_tot = 6.5 + nu_range = (5.0, np.inf) + reference = "10.1088/1361-6455/ab3c26" + + inner_channels = [ + AngularKetLS(l_c=0, l_r=1, l_tot=1, s_tot=1, j_tot=2, f_tot=6.5, species="Sr87"), # "5snp 3P2" + ] + outer_channels = [ + AngularKetFJ(l_c=0, l_r=1, j_c=0.5, f_c=5, j_r=1.5, f_tot=6.5, species="Sr87"), + ] + + eigen_quantum_defects = [ + [2.882, 0.446, -1.9], + ] + mixing_angles = [] + + +# -------------------------------------------------------- +# High-n D models +# -------------------------------------------------------- + + +class Sr87_D15_HighN(FModel): + species_name = "Sr87_mqdt" + name = "D F=3/2, nu > 25" + f_tot = 1.5 + nu_range = (25.0, np.inf) + reference = "10.1088/1361-6455/ab3c26" + + inner_channels = [ + AngularKetLS(l_c=0, l_r=2, l_tot=2, s_tot=1, j_tot=3, f_tot=1.5, species="Sr87"), # "5snd 3D3" + ] + outer_channels = [ + AngularKetFJ(l_c=0, l_r=2, j_c=0.5, f_c=4, j_r=2.5, f_tot=1.5, species="Sr87"), + ] + + eigen_quantum_defects = [ + [2.655, -41.4, -15363], + ] + mixing_angles = [] + + +class Sr87_D25_HighN(FModel): + species_name = "Sr87_mqdt" + name = "D F=5/2, nu > 25" + f_tot = 2.5 + nu_range = (25.0, np.inf) + reference = "10.1088/1361-6455/ab3c26" + + inner_channels = [ + AngularKetLS(l_c=0, l_r=2, l_tot=2, s_tot=0, j_tot=2, f_tot=2.5, species="Sr87"), # "5snd 1D2" + AngularKetLS(l_c=0, l_r=2, l_tot=2, s_tot=1, j_tot=2, f_tot=2.5, species="Sr87"), # "5snd 3D2" + AngularKetLS(l_c=0, l_r=2, l_tot=2, s_tot=1, j_tot=3, f_tot=2.5, species="Sr87"), # "5snd 3D3" + ] + outer_channels = [ + AngularKetFJ(l_c=0, l_r=2, j_c=0.5, f_c=4, j_r=1.5, f_tot=2.5, species="Sr87"), + AngularKetFJ(l_c=0, l_r=2, j_c=0.5, f_c=4, j_r=2.5, f_tot=2.5, species="Sr87"), + AngularKetFJ(l_c=0, l_r=2, j_c=0.5, f_c=5, j_r=2.5, f_tot=2.5, species="Sr87"), + ] + + eigen_quantum_defects = [ + [2.3847, -39.41, -1090], + [2.66149, -16.77, -6656], + [2.655, -41.4, -15363], + ] + mixing_angles = [ + (0, 1, -0.14), + ] + + +class Sr87_D35_HighN(FModel): + species_name = "Sr87_mqdt" + name = "D F=7/2, nu > 25" + f_tot = 3.5 + nu_range = (25.0, np.inf) + reference = "10.1088/1361-6455/ab3c26" + + inner_channels = [ + AngularKetLS(l_c=0, l_r=2, l_tot=2, s_tot=0, j_tot=2, f_tot=3.5, species="Sr87"), # "5snd 1D2" + AngularKetLS(l_c=0, l_r=2, l_tot=2, s_tot=1, j_tot=1, f_tot=3.5, species="Sr87"), # "5snd 3D1" + AngularKetLS(l_c=0, l_r=2, l_tot=2, s_tot=1, j_tot=2, f_tot=3.5, species="Sr87"), # "5snd 3D2" + AngularKetLS(l_c=0, l_r=2, l_tot=2, s_tot=1, j_tot=3, f_tot=3.5, species="Sr87"), # "5snd 3D3" + ] + outer_channels = [ + AngularKetFJ(l_c=0, l_r=2, j_c=0.5, f_c=4, j_r=1.5, f_tot=3.5, species="Sr87"), + AngularKetFJ(l_c=0, l_r=2, j_c=0.5, f_c=4, j_r=2.5, f_tot=3.5, species="Sr87"), + AngularKetFJ(l_c=0, l_r=2, j_c=0.5, f_c=5, j_r=1.5, f_tot=3.5, species="Sr87"), + AngularKetFJ(l_c=0, l_r=2, j_c=0.5, f_c=5, j_r=2.5, f_tot=3.5, species="Sr87"), + ] + + eigen_quantum_defects = [ + [2.3847, -39.41, -1090], + [2.67524, -13.15, -4444], + [2.66149, -16.77, -6656], + [2.655, -41.4, -15363], + ] + mixing_angles = [ + (0, 2, -0.14), + ] + + +class Sr87_D45_HighN(FModel): + species_name = "Sr87_mqdt" + name = "D F=9/2, nu > 25" + f_tot = 4.5 + nu_range = (25.0, np.inf) + reference = "10.1088/1361-6455/ab3c26" + + inner_channels = [ + AngularKetLS(l_c=0, l_r=2, l_tot=2, s_tot=0, j_tot=2, f_tot=4.5, species="Sr87"), # "5snd 1D2" + AngularKetLS(l_c=0, l_r=2, l_tot=2, s_tot=1, j_tot=1, f_tot=4.5, species="Sr87"), # "5snd 3D1" + AngularKetLS(l_c=0, l_r=2, l_tot=2, s_tot=1, j_tot=2, f_tot=4.5, species="Sr87"), # "5snd 3D2" + AngularKetLS(l_c=0, l_r=2, l_tot=2, s_tot=1, j_tot=3, f_tot=4.5, species="Sr87"), # "5snd 3D3" + ] + outer_channels = [ + AngularKetFJ(l_c=0, l_r=2, j_c=0.5, f_c=4, j_r=1.5, f_tot=4.5, species="Sr87"), + AngularKetFJ(l_c=0, l_r=2, j_c=0.5, f_c=4, j_r=2.5, f_tot=4.5, species="Sr87"), + AngularKetFJ(l_c=0, l_r=2, j_c=0.5, f_c=5, j_r=1.5, f_tot=4.5, species="Sr87"), + AngularKetFJ(l_c=0, l_r=2, j_c=0.5, f_c=5, j_r=2.5, f_tot=4.5, species="Sr87"), + ] + + eigen_quantum_defects = [ + [2.3847, -39.41, -1090], + [2.67524, -13.15, -4444], + [2.66149, -16.77, -6656], + [2.655, -41.4, -15363], + ] + mixing_angles = [ + (0, 2, -0.14), + ] + + +class Sr87_D55_HighN(FModel): + species_name = "Sr87_mqdt" + name = "D F=11/2, nu > 25" + f_tot = 5.5 + nu_range = (25.0, np.inf) + reference = "10.1088/1361-6455/ab3c26" + + inner_channels = [ + AngularKetLS(l_c=0, l_r=2, l_tot=2, s_tot=0, j_tot=2, f_tot=5.5, species="Sr87"), # "5snd 1D2" + AngularKetLS(l_c=0, l_r=2, l_tot=2, s_tot=1, j_tot=1, f_tot=5.5, species="Sr87"), # "5snd 3D1" + AngularKetLS(l_c=0, l_r=2, l_tot=2, s_tot=1, j_tot=2, f_tot=5.5, species="Sr87"), # "5snd 3D2" + AngularKetLS(l_c=0, l_r=2, l_tot=2, s_tot=1, j_tot=3, f_tot=5.5, species="Sr87"), # "5snd 3D3" + ] + outer_channels = [ + AngularKetFJ(l_c=0, l_r=2, j_c=0.5, f_c=4, j_r=1.5, f_tot=5.5, species="Sr87"), + AngularKetFJ(l_c=0, l_r=2, j_c=0.5, f_c=4, j_r=2.5, f_tot=5.5, species="Sr87"), + AngularKetFJ(l_c=0, l_r=2, j_c=0.5, f_c=5, j_r=1.5, f_tot=5.5, species="Sr87"), + AngularKetFJ(l_c=0, l_r=2, j_c=0.5, f_c=5, j_r=2.5, f_tot=5.5, species="Sr87"), + ] + + eigen_quantum_defects = [ + [2.3847, -39.41, -1090], + [2.67524, -13.15, -4444], + [2.66149, -16.77, -6656], + [2.655, -41.4, -15363], + ] + mixing_angles = [ + (0, 2, -0.14), + ] + + +class Sr87_D65_HighN(FModel): + species_name = "Sr87_mqdt" + name = "D F=13/2, nu > 25" + f_tot = 6.5 + nu_range = (25.0, np.inf) + reference = "10.1088/1361-6455/ab3c26" + + inner_channels = [ + AngularKetLS(l_c=0, l_r=2, l_tot=2, s_tot=0, j_tot=2, f_tot=6.5, species="Sr87"), # "5snd 1D2" + AngularKetLS(l_c=0, l_r=2, l_tot=2, s_tot=1, j_tot=2, f_tot=6.5, species="Sr87"), # "5snd 3D2" + AngularKetLS(l_c=0, l_r=2, l_tot=2, s_tot=1, j_tot=3, f_tot=6.5, species="Sr87"), # "5snd 3D3" + ] + outer_channels = [ + AngularKetFJ(l_c=0, l_r=2, j_c=0.5, f_c=4, j_r=2.5, f_tot=6.5, species="Sr87"), + AngularKetFJ(l_c=0, l_r=2, j_c=0.5, f_c=5, j_r=1.5, f_tot=6.5, species="Sr87"), + AngularKetFJ(l_c=0, l_r=2, j_c=0.5, f_c=5, j_r=2.5, f_tot=6.5, species="Sr87"), + ] + + eigen_quantum_defects = [ + [2.3847, -39.41, -1090], + [2.66149, -16.77, -6656], + [2.655, -41.4, -15363], + ] + mixing_angles = [ + (0, 1, -0.14), + ] + + +class Sr87_D75_HighN(FModel): + species_name = "Sr87_mqdt" + name = "D F=15/2, nu > 25" + f_tot = 7.5 + nu_range = (25.0, np.inf) + reference = "10.1088/1361-6455/ab3c26" + + inner_channels = [ + AngularKetLS(l_c=0, l_r=2, l_tot=2, s_tot=1, j_tot=3, f_tot=7.5, species="Sr87"), # "5snd 3D3" + ] + outer_channels = [ + AngularKetFJ(l_c=0, l_r=2, j_c=0.5, f_c=5, j_r=2.5, f_tot=7.5, species="Sr87"), + ] + + eigen_quantum_defects = [ + [2.655, -41.4, -15363], + ] + mixing_angles = [] + + +# -------------------------------------------------------- +# High-n F models +# -------------------------------------------------------- + + +class Sr87_F45_HighN(FModel): + species_name = "Sr87_mqdt" + name = "F F=9/2, nu > 9" + f_tot = 4.5 + nu_range = (9.0, np.inf) + reference = "10.1088/1361-6455/ab3c26" + + inner_channels = [ + AngularKetLS(l_c=0, l_r=3, l_tot=3, s_tot=0, j_tot=3, f_tot=4.5, species="Sr87"), # "5snf 1F3" + AngularKetLS(l_c=0, l_r=3, l_tot=3, s_tot=1, j_tot=2, f_tot=4.5, species="Sr87"), # "5snf 3F2" + AngularKetLS(l_c=0, l_r=3, l_tot=3, s_tot=1, j_tot=3, f_tot=4.5, species="Sr87"), # "5snf 3F3" + AngularKetLS(l_c=0, l_r=3, l_tot=3, s_tot=1, j_tot=4, f_tot=4.5, species="Sr87"), # "5snf 3F4" + ] + outer_channels = [ + AngularKetFJ(l_c=0, l_r=3, j_c=0.5, f_c=4, j_r=2.5, f_tot=4.5, species="Sr87"), + AngularKetFJ(l_c=0, l_r=3, j_c=0.5, f_c=4, j_r=3.5, f_tot=4.5, species="Sr87"), + AngularKetFJ(l_c=0, l_r=3, j_c=0.5, f_c=5, j_r=2.5, f_tot=4.5, species="Sr87"), + AngularKetFJ(l_c=0, l_r=3, j_c=0.5, f_c=5, j_r=3.5, f_tot=4.5, species="Sr87"), + ] + + eigen_quantum_defects = [ + [0.12, -2.2, 120], + [0.12, -2.2, 120], + [0.12, -2.2, 120], + [0.12, -2.2, 120], + ] + mixing_angles = [] diff --git a/src/rydstate/species/mqdt/sr88.py b/src/rydstate/species/mqdt/sr88.py new file mode 100644 index 0000000..18ed1bc --- /dev/null +++ b/src/rydstate/species/mqdt/sr88.py @@ -0,0 +1,303 @@ +from __future__ import annotations + +import numpy as np + +from rydstate.angular.angular_ket import AngularKetFJ, AngularKetLS +from rydstate.angular.core_ket_base import CoreKet +from rydstate.species.mqdt.fmodel import FModel +from rydstate.species.mqdt.species_object_mqdt import SpeciesObjectMQDT +from rydstate.units import electron_mass, rydberg_constant + + +class Strontium88MQDT(SpeciesObjectMQDT): + name = "Sr88_mqdt" + Z = 38 + i_c = 0 + number_valence_electrons = 2 + + potential_type_default = "model_potential_fei_2009" + + # https://physics.nist.gov/PhysRefData/Handbook/Tables/strontiumtable1.htm + _isotope_mass = 87.9056122571 # u + _corrected_rydberg_constant = ( + rydberg_constant.m / (1 + electron_mass.to("u").m / _isotope_mass), + None, + str(rydberg_constant.u), + ) + + _ionization_threshold_dict = { + CoreKet(i_c, 0.5, 0, 0.5): (45932.1956, None, "1/cm"), + } + core_ground_state = CoreKet(i_c, 0.5, 0, 0.5) + nuclear_dipole = 2.3 + + +# -------------------------------------------------------- +# MQDT models valid at large n +# -------------------------------------------------------- + + +class Sr88_S0_HighN(FModel): + species_name = "Sr88_mqdt" + name = "S J=0, nu > 10" + f_tot = 0 + nu_range = (10.0, np.inf) + reference = "10.1088/1361-6455/ab3c26" + + inner_channels = [ + AngularKetLS(l_c=0, l_r=0, l_tot=0, s_tot=0, j_tot=0, species="Sr88"), # "5sns 1S0" + ] + outer_channels = [ + AngularKetFJ(l_c=0, l_r=0, j_c=0.5, j_r=0.5, f_tot=0, species="Sr88"), + ] + + eigen_quantum_defects = [ + [3.26896, -0.138, 0.9], + ] + mixing_angles = [] + + +class Sr88_S1_HighN(FModel): + species_name = "Sr88_mqdt" + name = "S J=1, nu > 11" + f_tot = 1 + nu_range = (11.0, np.inf) + reference = "10.1088/1361-6455/ab3c26" + + inner_channels = [ + AngularKetLS(l_c=0, l_r=0, l_tot=0, s_tot=1, j_tot=1, species="Sr88"), # "5sns 3S1" + ] + outer_channels = [ + AngularKetFJ(l_c=0, l_r=0, j_c=0.5, j_r=0.5, f_tot=1, species="Sr88"), + ] + + eigen_quantum_defects = [ + [3.370778, 0.418, -0.3], + ] + mixing_angles = [] + + +# -------------------------------------------------------- +# MQDT models valid at small n +# -------------------------------------------------------- + + +class Sr88_P1_LowN(FModel): + species_name = "Sr88_mqdt" + name = "P J=1 (recombination), 1.8 < nu < 2.2" + f_tot = 1 + nu_range = (1.8, 2.2) + reference = "10.1088/1361-6455/ab3c26" + + inner_channels = [ + AngularKetLS(l_c=0, l_r=1, l_tot=1, s_tot=0, j_tot=1, species="Sr88"), # "5snp 1P1" + AngularKetLS(l_c=0, l_r=1, l_tot=1, s_tot=1, j_tot=1, species="Sr88"), # "5snp 3P1" + ] + outer_channels = [ + AngularKetFJ(l_c=0, l_r=1, j_c=0.5, j_r=0.5, f_tot=1, species="Sr88"), + AngularKetFJ(l_c=0, l_r=1, j_c=0.5, j_r=1.5, f_tot=1, species="Sr88"), + ] + + eigen_quantum_defects = [ + [0.87199081, 0], + [0.13140955, 0], + ] + mixing_angles = [ + (0, 1, [1.31169947, -4.48280597]), + ] + + +# -------------------------------------------------------- +# MQDT models valid at large n (continued) +# -------------------------------------------------------- + + +class Sr88_P0_HighN(FModel): + species_name = "Sr88_mqdt" + name = "P J=0, nu > 7" + f_tot = 0 + nu_range = (7.0, np.inf) + reference = "10.1088/1361-6455/ab3c26" + + inner_channels = [ + AngularKetLS(l_c=0, l_r=1, l_tot=1, s_tot=1, j_tot=0, species="Sr88"), # "5snp 3P0" + ] + outer_channels = [ + AngularKetFJ(l_c=0, l_r=1, j_c=0.5, j_r=0.5, f_tot=0, species="Sr88"), + ] + + eigen_quantum_defects = [ + [2.8867, 0.44, -1.9], + ] + mixing_angles = [] + + +class Sr88_P1_HighN(FModel): + species_name = "Sr88_mqdt" + name = "P J=1, nu > 5" + f_tot = 1 + nu_range = (5.0, np.inf) + reference = "10.1088/1361-6455/ab3c26" + + inner_channels = [ + AngularKetLS(l_c=0, l_r=1, l_tot=1, s_tot=0, j_tot=1, species="Sr88"), # "5snp 1P1" + AngularKetLS(l_c=0, l_r=1, l_tot=1, s_tot=1, j_tot=1, species="Sr88"), # "5snp 3P1" + ] + outer_channels = [ + AngularKetFJ(l_c=0, l_r=1, j_c=0.5, j_r=0.5, f_tot=1, species="Sr88"), + AngularKetFJ(l_c=0, l_r=1, j_c=0.5, j_r=1.5, f_tot=1, species="Sr88"), + ] + + eigen_quantum_defects = [ + [2.724, -4.67, -157], + [2.8826, 0.407, -1.3], + ] + mixing_angles = [] + + +class Sr88_P2_HighN(FModel): + species_name = "Sr88_mqdt" + name = "P J=2, nu > 5" + f_tot = 2 + nu_range = (5.0, np.inf) + reference = "10.1088/1361-6455/ab3c26" + + inner_channels = [ + AngularKetLS(l_c=0, l_r=1, l_tot=1, s_tot=1, j_tot=2, species="Sr88"), # "5snp 3P2" + ] + outer_channels = [ + AngularKetFJ(l_c=0, l_r=1, j_c=0.5, j_r=1.5, f_tot=2, species="Sr88"), + ] + + eigen_quantum_defects = [ + [2.882, 0.446, -1.9], + ] + mixing_angles = [] + + +class Sr88_D1_HighN(FModel): + species_name = "Sr88_mqdt" + name = "D J=1, nu > 17" + f_tot = 1 + nu_range = (17.0, np.inf) + reference = "10.1088/1361-6455/ab3c26" + + inner_channels = [ + AngularKetLS(l_c=0, l_r=2, l_tot=2, s_tot=1, j_tot=1, species="Sr88"), # "5snd 3D1" + ] + outer_channels = [ + AngularKetFJ(l_c=0, l_r=2, j_c=0.5, j_r=1.5, f_tot=1, species="Sr88"), + ] + + eigen_quantum_defects = [ + [2.67524, -13.15, -4444], + ] + mixing_angles = [] + + +class Sr88_D2_HighN(FModel): + species_name = "Sr88_mqdt" + name = "D J=2, nu > 25" + f_tot = 2 + nu_range = (25.0, np.inf) + reference = "10.1088/1361-6455/ab3c26" + + inner_channels = [ + AngularKetLS(l_c=0, l_r=2, l_tot=2, s_tot=0, j_tot=2, species="Sr88"), # "5snd 1D2" + AngularKetLS(l_c=0, l_r=2, l_tot=2, s_tot=1, j_tot=2, species="Sr88"), # "5snd 3D2" + ] + outer_channels = [ + AngularKetFJ(l_c=0, l_r=2, j_c=0.5, j_r=1.5, f_tot=2, species="Sr88"), + AngularKetFJ(l_c=0, l_r=2, j_c=0.5, j_r=2.5, f_tot=2, species="Sr88"), + ] + + eigen_quantum_defects = [ + [2.3847, -39.41, -1090], + [2.66149, -16.77, -6656], + ] + mixing_angles = [ + (0, 1, -0.14), + ] + + +class Sr88_D3_HighN(FModel): + species_name = "Sr88_mqdt" + name = "D J=3, nu > 25" + f_tot = 3 + nu_range = (25.0, np.inf) + reference = "10.1088/1361-6455/ab3c26" + + inner_channels = [ + AngularKetLS(l_c=0, l_r=2, l_tot=2, s_tot=1, j_tot=3, species="Sr88"), # "5snd 3D3" + ] + outer_channels = [ + AngularKetFJ(l_c=0, l_r=2, j_c=0.5, j_r=2.5, f_tot=3, species="Sr88"), + ] + + eigen_quantum_defects = [ + [2.655, -41.4, -15363], + ] + mixing_angles = [] + + +class Sr88_F2_HighN(FModel): + species_name = "Sr88_mqdt" + name = "F J=2, nu > 9" + f_tot = 2 + nu_range = (9.0, np.inf) + reference = "10.1088/1361-6455/ab3c26" + + inner_channels = [ + AngularKetLS(l_c=0, l_r=3, l_tot=3, s_tot=1, j_tot=2, species="Sr88"), # "5snf 3F2" + ] + outer_channels = [ + AngularKetFJ(l_c=0, l_r=3, j_c=0.5, j_r=2.5, f_tot=2, species="Sr88"), + ] + + eigen_quantum_defects = [ + [0.12, -2.2, 120], + ] + mixing_angles = [] + + +class Sr88_F3_HighN(FModel): + species_name = "Sr88_mqdt" + name = "F J=3, nu > 9" + f_tot = 3 + nu_range = (9.0, np.inf) + reference = "10.1088/1361-6455/ab3c26" + + inner_channels = [ + AngularKetLS(l_c=0, l_r=3, l_tot=3, s_tot=0, j_tot=3, species="Sr88"), # "5snf 1F3" + AngularKetLS(l_c=0, l_r=3, l_tot=3, s_tot=1, j_tot=3, species="Sr88"), # "5snf 3F3" + ] + outer_channels = [ + AngularKetFJ(l_c=0, l_r=3, j_c=0.5, j_r=2.5, f_tot=3, species="Sr88"), + AngularKetFJ(l_c=0, l_r=3, j_c=0.5, j_r=3.5, f_tot=3, species="Sr88"), + ] + + eigen_quantum_defects = [ + [0.089, -2, 30], + [0.12, -2.2, 120], + ] + mixing_angles = [] + + +class Sr88_F4_HighN(FModel): + species_name = "Sr88_mqdt" + name = "F J=4, nu > 9" + f_tot = 4 + nu_range = (9.0, np.inf) + reference = "10.1088/1361-6455/ab3c26" + + inner_channels = [ + AngularKetLS(l_c=0, l_r=3, l_tot=3, s_tot=1, j_tot=4, species="Sr88"), # "5snf 3F4" + ] + outer_channels = [ + AngularKetFJ(l_c=0, l_r=3, j_c=0.5, j_r=3.5, f_tot=4, species="Sr88"), + ] + + eigen_quantum_defects = [ + [0.12, -2.2, 120], + ] + mixing_angles = [] diff --git a/src/rydstate/species/mqdt/yb171.py b/src/rydstate/species/mqdt/yb171.py new file mode 100644 index 0000000..4de7fdf --- /dev/null +++ b/src/rydstate/species/mqdt/yb171.py @@ -0,0 +1,1043 @@ +from __future__ import annotations + +import numpy as np + +from rydstate.angular.angular_ket import AngularKetFJ, AngularKetJJ, AngularKetLS +from rydstate.angular.angular_ket_dummy import AngularKetDummy +from rydstate.angular.core_ket_base import CoreKet, CoreKetDummy +from rydstate.angular.utils import Unknown +from rydstate.species.mqdt.fmodel import FModel +from rydstate.species.mqdt.species_object_mqdt import SpeciesObjectMQDT +from rydstate.units import electron_mass, rydberg_constant + + +class Ytterbium171MQDT(SpeciesObjectMQDT): + name = "Yb171_mqdt" + Z = 70 + i_c = 0.5 + number_valence_electrons = 2 + + potential_type_default = "model_potential_fei_2009" + + # https://iopscience.iop.org/article/10.1088/1674-1056/18/10/025 + model_potential_parameter_fei_2009 = (0.8704, 22.0040, 0.1513, 0.3306) + + # https://physics.nist.gov/PhysRefData/Handbook/Tables/ytterbiumtable1.htm + _isotope_mass = 170.9363258 # u + _corrected_rydberg_constant = ( + rydberg_constant.m / (1 + electron_mass.to("u").m / _isotope_mass), + None, + str(rydberg_constant.u), + ) + + _ionization_threshold_dict = { + CoreKet(i_c, 0.5, 0, 0.5, 0): (50442.795744, None, "1/cm"), + CoreKet(i_c, 0.5, 0, 0.5, 1): (50443.217463, None, "1/cm"), + CoreKet(i_c, 0.5, 1, 0.5, Unknown): (77504.98, None, "1/cm"), + CoreKet(i_c, 0.5, 1, Unknown, Unknown): (79725.35, None, "1/cm"), + CoreKet(i_c, 0.5, 1, 1.5, Unknown): (80835.39, None, "1/cm"), + CoreKetDummy("4f13 5d 6s"): (83967.7, None, "1/cm"), + } + core_ground_state = CoreKet(i_c, 0.5, 0, 0.5, 1) + nuclear_dipole = 0.49367 + + +# -------------------------------------------------------- +# MQDT models valid at large n +# -------------------------------------------------------- + + +class Yb171_S05_HighN(FModel): + species_name = "Yb171_mqdt" + name = "S F=0.5, nu > 26" + f_tot = 0.5 + nu_range = (26.0, np.inf) + reference = "10.1103/PhysRevA.112.042817, 10.1103/PhysRevX.15.011009" + + inner_channels = [ + AngularKetLS(l_c=0, l_r=0, l_tot=0, s_tot=0, j_tot=0, f_tot=0.5, species="Yb171"), # "6sns 1S0" + AngularKetDummy("4f13 5d 6snl a", f_tot=0.5), + AngularKetLS(l_c=1, l_r=1, l_tot=0, s_tot=0, j_tot=0, f_tot=0.5, species="Yb171"), # "6pnp 1S0" + AngularKetDummy("4f13 5d 6snl b", f_tot=0.5), + AngularKetLS(l_c=1, l_r=1, l_tot=1, s_tot=1, j_tot=0, f_tot=0.5, species="Yb171"), # "6pnp 3P0" + AngularKetDummy("4f13 5d 6snl c", f_tot=0.5), + AngularKetLS(l_c=0, l_r=0, l_tot=0, s_tot=1, j_tot=1, f_tot=0.5, species="Yb171"), # "6sns 3S1" + ] + outer_channels = [ + AngularKetFJ(l_c=0, l_r=0, j_c=0.5, f_c=0, j_r=0.5, f_tot=0.5, species="Yb171"), + AngularKetDummy("4f13 5d 6snl a", f_tot=0.5), + AngularKetFJ(l_c=1, l_r=1, j_c=1.5, f_c=1, j_r=1.5, f_tot=0.5, species="Yb171"), + AngularKetDummy("4f13 5d 6snl b", f_tot=0.5), + AngularKetFJ(l_c=1, l_r=1, j_c=0.5, f_c=Unknown, j_r=0.5, f_tot=0.5, species="Yb171"), + AngularKetDummy("4f13 5d 6snl c", f_tot=0.5), + AngularKetFJ(l_c=0, l_r=0, j_c=0.5, f_c=1, j_r=0.5, f_tot=0.5, species="Yb171"), + ] + + eigen_quantum_defects = [ + [0.357488757, 0.165981371, 0, 0, 0], + [0.203918644, 0, 0, 0, 0], + [0.116819032, 0, 0, 0, 0], + [0.287350241, 0, 0, 0, 0], + [0.247621114, 0, 0, 0, 0], + [0.148681324, 0, 0, 0, 0], + [0.438542187, 3.78366407, -10709.7378, 8054542.58, -2523011670], + ] + mixing_angles = [ + (0, 1, 0.131755467), + (0, 2, 0.297504211), + (0, 3, 0.055421439), + (2, 3, 0.100871756), + (2, 4, 0.103123032), + (0, 5, 0.137753117), + ] + manual_frame_transformation_outer_inner = np.array( + [ + [1 / 2, 0, 0, 0, 0, 0, np.sqrt(3) / 2], + [0, 1, 0, 0, 0, 0, 0], + [0, 0, np.sqrt(2 / 3), 0, -np.sqrt(1 / 3), 0, 0], + [0, 0, 0, 1, 0, 0, 0], + [0, 0, np.sqrt(1 / 3), 0, np.sqrt(2 / 3), 0, 0], + [0, 0, 0, 0, 0, 1, 0], + [np.sqrt(3) / 2, 0, 0, 0, 0, 0, -1 / 2], + ] + ) + + +class Yb171_S15_HighN(FModel): + species_name = "Yb171_mqdt" + name = "S F=1.5, nu > 26" + f_tot = 1.5 + nu_range = (26.0, np.inf) + reference = "10.1103/PhysRevX.15.011009" + + inner_channels = [ + AngularKetLS(l_c=0, l_r=0, l_tot=0, s_tot=1, j_tot=1, f_tot=1.5, species="Yb171"), # "6sns 3S1" + ] + outer_channels = [ + AngularKetFJ(l_c=0, l_r=0, j_c=0.5, f_c=1, j_r=0.5, f_tot=1.5, species="Yb171"), + ] + + eigen_quantum_defects = [ + [0.438426851, 3.91762642, -10612.6828, 8017432.38, -2582622910.0], + ] + mixing_angles = [] + + +class Yb171_P05_HighN(FModel): + species_name = "Yb171_mqdt" + name = "P F=0.5, nu > 5.7" + f_tot = 0.5 + nu_range = (5.7, np.inf) + reference = "10.1103/PhysRevA.112.042817, 10.1103/PhysRevX.15.011009" + + inner_channels = [ + AngularKetLS(l_c=0, l_r=1, l_tot=1, s_tot=0, j_tot=1, f_tot=0.5, species="Yb171"), # "6snp 1P1" + AngularKetLS(l_c=0, l_r=1, l_tot=1, s_tot=1, j_tot=1, f_tot=0.5, species="Yb171"), # "6snp 3P1" + AngularKetDummy("4f13 5d 6snl a", f_tot=0.5), + AngularKetDummy("4f13 5d 6snl b", f_tot=0.5), + AngularKetDummy("4f13 5d 6snl c", f_tot=0.5), + AngularKetDummy("4f13 5d 6snl d", f_tot=0.5), + AngularKetLS(l_c=0, l_r=1, l_tot=1, s_tot=1, j_tot=0, f_tot=0.5, species="Yb171"), # "6snp 3P0" + AngularKetDummy("4f13 5d 6snl e", f_tot=0.5), + ] + outer_channels = [ + AngularKetFJ(l_c=0, l_r=1, j_c=0.5, f_c=1, j_r=1.5, f_tot=0.5, species="Yb171"), + AngularKetFJ(l_c=0, l_r=1, j_c=0.5, f_c=1, j_r=0.5, f_tot=0.5, species="Yb171"), + AngularKetDummy("4f13 5d 6snl a", f_tot=0.5), + AngularKetDummy("4f13 5d 6snl b", f_tot=0.5), + AngularKetDummy("4f13 5d 6snl c", f_tot=0.5), + AngularKetDummy("4f13 5d 6snl d", f_tot=0.5), + AngularKetFJ(l_c=0, l_r=1, j_c=0.5, f_c=0, j_r=0.5, f_tot=0.5, species="Yb171"), + AngularKetDummy("4f13 5d 6snl e", f_tot=0.5), + ] + + eigen_quantum_defects = [ + [0.922094502, 2.12370136], + [0.981191543, -4.54209175], + [0.229094016, 0], + [0.206073107, 0], + [0.193527627, 0], + [0.181165673, 0], + [0.953185132, 0.0277444042], + [0.198448494, 0], + ] + mixing_angles = [ + (0, 1, [-0.102285383, 153.521338, -15393.2283]), + (1, 6, -0.00168607392), + (0, 2, -0.0719467433), + (0, 3, -0.0673315968), + (0, 4, -0.0221077377), + (0, 5, -0.107638329), + (1, 2, 0.0416653549), + (1, 3, 0.0590660991), + (1, 4, 0.0861585559), + (1, 5, 0.0566417469), + (6, 7, 0.163113423), + ] + + +class Yb171_P15_HighN(FModel): + species_name = "Yb171_mqdt" + name = "P F=1.5, nu > 10" + f_tot = 1.5 + nu_range = (10.0, np.inf) + reference = "10.1103/PhysRevA.112.042817" + + inner_channels = [ + AngularKetLS(l_c=0, l_r=1, l_tot=1, s_tot=0, j_tot=1, f_tot=1.5, species="Yb171"), # "6snp 1P1" + AngularKetLS(l_c=0, l_r=1, l_tot=1, s_tot=1, j_tot=1, f_tot=1.5, species="Yb171"), # "6snp 3P1" + AngularKetDummy("4f13 5d 6snl a", f_tot=1.5), + AngularKetDummy("4f13 5d 6snl b", f_tot=1.5), + AngularKetDummy("4f13 5d 6snl c", f_tot=1.5), + AngularKetDummy("4f13 5d 6snl d", f_tot=1.5), + AngularKetLS(l_c=0, l_r=1, l_tot=1, s_tot=1, j_tot=2, f_tot=1.5, species="Yb171"), # "6snp 3P2" + AngularKetDummy("4f13 5d 6snl e", f_tot=1.5), + AngularKetDummy("4f13 5d 6snl f", f_tot=1.5), + AngularKetDummy("4f13 5d 6snl g", f_tot=1.5), + AngularKetLS(l_c=0, l_r=3, l_tot=3, s_tot=1, j_tot=2, f_tot=1.5, species="Yb171"), # "6snf 3F2" + ] + outer_channels = [ + AngularKetFJ(l_c=0, l_r=1, j_c=0.5, f_c=1, j_r=1.5, f_tot=1.5, species="Yb171"), + AngularKetFJ(l_c=0, l_r=1, j_c=0.5, f_c=1, j_r=0.5, f_tot=1.5, species="Yb171"), + AngularKetDummy("4f13 5d 6snl a", f_tot=1.5), + AngularKetDummy("4f13 5d 6snl b", f_tot=1.5), + AngularKetDummy("4f13 5d 6snl c", f_tot=1.5), + AngularKetDummy("4f13 5d 6snl d", f_tot=1.5), + AngularKetFJ(l_c=0, l_r=1, j_c=0.5, f_c=0, j_r=1.5, f_tot=1.5, species="Yb171"), + AngularKetDummy("4f13 5d 6snl e", f_tot=1.5), + AngularKetDummy("4f13 5d 6snl f", f_tot=1.5), + AngularKetDummy("4f13 5d 6snl g", f_tot=1.5), + AngularKetFJ(l_c=0, l_r=3, j_c=0.5, f_c=1, j_r=2.5, f_tot=1.5, species="Yb171"), + ] + + eigen_quantum_defects = [ + [0.922094502, 2.12370136, 0], + [0.981191543, -4.54209175, 0], + [0.229094016, 0, 0], + [0.206073107, 0, 0], + [0.193527627, 0, 0], + [0.181165673, 0, 0], + [0.925345494, -3.23594086, 80.2535181], + [0.232649227, 0, 0], + [0.210070444, 0, 0], + [0.185699031, 0, 0], + [0.0718955585, -1.0913707, -38.4618954], + ] + mixing_angles = [ + (0, 1, [-0.102285383, 153.251338, -15393.2283]), + (0, 2, -0.0719467433), + (0, 3, -0.0673315968), + (0, 4, -0.0221077377), + (0, 5, -0.107638329), + (1, 2, 0.0416653549), + (1, 3, 0.0590660991), + (1, 4, 0.0861585559), + (1, 5, 0.0566417469), + (6, 7, 0.0703574701), + (6, 8, 0.0235308506), + (6, 9, -0.0295876723), + (6, 10, 0.018377516), + ] + + +class Yb171_D05_HighN(FModel): + species_name = "Yb171_mqdt" + name = "D F=0.5, nu > 30" + f_tot = 0.5 + nu_range = (30.0, np.inf) + reference = "10.1103/PhysRevX.15.011009" + + inner_channels = [ + AngularKetLS(l_c=0, l_r=2, l_tot=2, s_tot=1, j_tot=1, f_tot=0.5, species="Yb171"), # "6snd 3D1" + ] + outer_channels = [ + AngularKetFJ(l_c=0, l_r=2, j_c=0.5, f_c=1, j_r=1.5, f_tot=0.5, species="Yb171"), + ] + + eigen_quantum_defects = [ + [0.75258093, 0.382628525, -483.120633], + ] + mixing_angles = [] + + +class Yb171_D15_HighN(FModel): + species_name = "Yb171_mqdt" + name = "D F=1.5, nu > 30" + f_tot = 1.5 + nu_range = (30.0, np.inf) + reference = "10.1103/PhysRevA.112.042817, 10.1103/PhysRevX.15.011009" + + inner_channels = [ + AngularKetLS(l_c=0, l_r=2, l_tot=2, s_tot=0, j_tot=2, f_tot=1.5, species="Yb171"), # "6snd 1D2" + AngularKetLS(l_c=0, l_r=2, l_tot=2, s_tot=1, j_tot=2, f_tot=1.5, species="Yb171"), # "6snd 3D2" + AngularKetDummy("4f13 5d 6snl a", f_tot=1.5), + AngularKetDummy("4f13 5d 6snl b", f_tot=1.5), + AngularKetLS(l_c=1, l_r=1, l_tot=2, s_tot=0, j_tot=2, f_tot=1.5, species="Yb171"), # "6pnp 1D2" + AngularKetLS(l_c=0, l_r=2, l_tot=2, s_tot=1, j_tot=1, f_tot=1.5, species="Yb171"), # "6snd 3D1" + ] + outer_channels = [ + AngularKetFJ(l_c=0, l_r=2, j_c=0.5, f_c=1, j_r=2.5, f_tot=1.5, species="Yb171"), + AngularKetFJ(l_c=0, l_r=2, j_c=0.5, f_c=1, j_r=1.5, f_tot=1.5, species="Yb171"), + AngularKetDummy("4f13 5d 6snl a", f_tot=1.5), + AngularKetDummy("4f13 5d 6snl b", f_tot=1.5), + AngularKetFJ(l_c=1, l_r=1, j_c=Unknown, f_c=Unknown, j_r=Unknown, f_tot=1.5, species="Yb171"), + AngularKetFJ(l_c=0, l_r=2, j_c=0.5, f_c=0, j_r=1.5, f_tot=1.5, species="Yb171"), + ] + manual_frame_transformation_outer_inner = np.array( + [ + [-np.sqrt(3 / 5), -np.sqrt(2 / 5), 0, 0, 0, 0], + [np.sqrt(3 / 5) / 2, -3 / (2 * np.sqrt(10)), 0, 0, 0, np.sqrt(5 / 2) / 2], + [0, 0, 1, 0, 0, 0], + [0, 0, 0, 1, 0, 0], + [0, 0, 0, 0, 1, 0], + [-1 / 2, np.sqrt(3 / 2) / 2, 0, 0, 0, np.sqrt(3 / 2) / 2], + ] + ) + + eigen_quantum_defects = [ + [0.73056016, -0.108286264, 0], + [0.75155852, 0.000367204397, 0], + [0.195831577, 0, 0], + [0.236133225, 0, 0], + [0.147506921, 0, 0], + [0.75336354, -1.84349555, 994.210321], + ] + mixing_angles = [ + (0, 1, [0.22146327, -16.2798928]), + (0, 2, 0.00431695191), + (0, 3, 0.0381576181), + (1, 3, -0.00708200703), + (0, 4, 0.109346659), + (1, 4, 0.0636016813), + ] + + +class Yb171_D25_HighN(FModel): + species_name = "Yb171_mqdt" + name = "D F=2.5, nu > 30" + f_tot = 2.5 + nu_range = (30.0, np.inf) + reference = "10.1103/PhysRevA.112.042817, 10.1103/PhysRevX.15.011009" + + inner_channels = [ + AngularKetLS(l_c=0, l_r=2, l_tot=2, s_tot=0, j_tot=2, f_tot=2.5, species="Yb171"), # "6snd 1D2" + AngularKetLS(l_c=0, l_r=2, l_tot=2, s_tot=1, j_tot=2, f_tot=2.5, species="Yb171"), # "6snd 3D2" + AngularKetDummy("4f13 5d 6snl a", f_tot=2.5), + AngularKetDummy("4f13 5d 6snl b", f_tot=2.5), + AngularKetLS(l_c=1, l_r=1, l_tot=2, s_tot=0, j_tot=2, f_tot=2.5, species="Yb171"), # "6pnp 1D2" + AngularKetLS(l_c=0, l_r=2, l_tot=2, s_tot=1, j_tot=3, f_tot=2.5, species="Yb171"), # "6snd 3D3" + ] + outer_channels = [ + AngularKetFJ(l_c=0, l_r=2, j_c=0.5, f_c=1, j_r=2.5, f_tot=2.5, species="Yb171"), + AngularKetFJ(l_c=0, l_r=2, j_c=0.5, f_c=1, j_r=1.5, f_tot=2.5, species="Yb171"), + AngularKetDummy("4f13 5d 6snl a", f_tot=2.5), + AngularKetDummy("4f13 5d 6snl b", f_tot=2.5), + AngularKetFJ(l_c=1, l_r=1, j_c=Unknown, f_c=Unknown, j_r=Unknown, f_tot=2.5, species="Yb171"), + AngularKetFJ(l_c=0, l_r=2, j_c=0.5, f_c=0, j_r=2.5, f_tot=2.5, species="Yb171"), + ] + manual_frame_transformation_outer_inner = np.array( + [ + [np.sqrt(7 / 5) / 2, np.sqrt(7 / 30), 0, 0, 0, -np.sqrt(5 / 3) / 2], + [-np.sqrt(2 / 5), np.sqrt(3 / 5), 0, 0, 0, 0], + [0, 0, 1, 0, 0, 0], + [0, 0, 0, 1, 0, 0], + [0, 0, 0, 0, 1, 0], + [1 / 2, np.sqrt(1 / 6), 0, 0, 0, np.sqrt(7 / 3) / 2], + ] + ) + + eigen_quantum_defects = [ + [0.73056016, -0.108286264, 0], + [0.75155852, 0.000367204397, 0], + [0.195831577, 0, 0], + [0.236133225, 0, 0], + [0.147506921, 0, 0], + [0.72861481, 0.79979111, -484.236631], + ] + mixing_angles = [ + (0, 1, [0.22146327, -16.2798928]), + (0, 2, 0.00431695191), + (0, 3, 0.0381576181), + (1, 3, -0.00708200703), + (0, 4, 0.109346659), + (1, 4, 0.0636016813), + ] + + +class Yb171_D35_HighN(FModel): + species_name = "Yb171_mqdt" + name = "D F=3.5, nu > 14" + f_tot = 3.5 + nu_range = (14.0, np.inf) + reference = "10.1103/PhysRevX.15.011009" + + inner_channels = [ + AngularKetLS(l_c=0, l_r=2, l_tot=2, s_tot=1, j_tot=3, f_tot=3.5, species="Yb171"), # "6snd 3D3" + ] + outer_channels = [ + AngularKetFJ(l_c=0, l_r=2, j_c=0.5, f_c=1, j_r=2.5, f_tot=3.5, species="Yb171"), + ] + + eigen_quantum_defects = [ + [0.72895315, -0.20653489, 220.484722], + ] + mixing_angles = [] + + +class Yb171_F25_HighN(FModel): + species_name = "Yb171_mqdt" + name = "F F=2.5, nu > 20" + f_tot = 2.5 + nu_range = (20.0, np.inf) + reference = "10.1103/PhysRevA.112.042817" + + inner_channels = [ + AngularKetLS(l_c=0, l_r=3, l_tot=3, s_tot=0, j_tot=3, f_tot=2.5, species="Yb171"), # "6snf 1F3" + AngularKetLS(l_c=0, l_r=3, l_tot=3, s_tot=1, j_tot=3, f_tot=2.5, species="Yb171"), # "6snf 3F3" + AngularKetDummy("4f13 5d 6snl a", f_tot=2.5), + AngularKetDummy("4f13 5d 6snl b", f_tot=2.5), + AngularKetDummy("4f13 5d 6snl c", f_tot=2.5), + AngularKetDummy("4f13 5d 6snl d", f_tot=2.5), + AngularKetDummy("4f13 5d 6snl e", f_tot=2.5), + AngularKetLS(l_c=0, l_r=1, l_tot=1, s_tot=1, j_tot=2, f_tot=2.5, species="Yb171"), # "6snf 3P2" + AngularKetDummy("4f13 5d 6snl f", f_tot=2.5), + AngularKetDummy("4f13 5d 6snl g", f_tot=2.5), + AngularKetDummy("4f13 5d 6snl h", f_tot=2.5), + AngularKetLS(l_c=0, l_r=3, l_tot=3, s_tot=1, j_tot=2, f_tot=2.5, species="Yb171"), # "6snf 3F2" + ] + outer_channels = [ + AngularKetFJ(l_c=0, l_r=3, j_c=0.5, f_c=1, j_r=3.5, f_tot=2.5, species="Yb171"), + AngularKetFJ(l_c=0, l_r=3, j_c=0.5, f_c=1, j_r=2.5, f_tot=2.5, species="Yb171"), + AngularKetDummy("4f13 5d 6snl a", f_tot=2.5), + AngularKetDummy("4f13 5d 6snl b", f_tot=2.5), + AngularKetDummy("4f13 5d 6snl c", f_tot=2.5), + AngularKetDummy("4f13 5d 6snl d", f_tot=2.5), + AngularKetDummy("4f13 5d 6snl e", f_tot=2.5), + AngularKetFJ(l_c=0, l_r=1, j_c=0.5, f_c=1, j_r=1.5, f_tot=2.5, species="Yb171"), + AngularKetDummy("4f13 5d 6snl f", f_tot=2.5), + AngularKetDummy("4f13 5d 6snl g", f_tot=2.5), + AngularKetDummy("4f13 5d 6snl h", f_tot=2.5), + AngularKetFJ(l_c=0, l_r=3, j_c=0.5, f_c=0, j_r=2.5, f_tot=2.5, species="Yb171"), + ] + + eigen_quantum_defects = [ + [0.277086649, -13.2829133, 0], + [0.0719837014, -0.741492064, 0], + [0.251457795, 0, 0], + [0.227434828, 0, 0], + [0.175780645, 0, 0], + [0.196547521, 0, 0], + [0.21440857, 0, 0], + [0.925345494, -3.23594086, 80.2535181], + [0.232649227, 0, 0], + [0.210070444, 0, 0], + [0.185699031, 0, 0], + [0.0718955585, -1.0913707, -38.4618954], + ] + mixing_angles = [ + (0, 1, [-0.0209955122, 0.251041249]), + (0, 2, -0.0585753224), + (0, 3, -0.0750574327), + (0, 4, 0.122671919), + (0, 5, -0.0401036164), + (0, 6, 0.0654271994), + (1, 2, -0.0683007974), + (1, 3, 0.035415976), + (1, 4, -0.0327625807), + (1, 5, -0.050225071), + (1, 6, 0.0455759316), + (7, 8, 0.0703574701), + (7, 9, 0.0235308506), + (7, 10, -0.0295876723), + (7, 11, 0.018377516), + ] + + +class Yb171_F35_HighN(FModel): + species_name = "Yb171_mqdt" + name = "F F=3.5, nu > 20" + f_tot = 3.5 + nu_range = (20.0, np.inf) + reference = "10.1103/PhysRevA.112.042817" + + inner_channels = [ + AngularKetLS(l_c=0, l_r=3, l_tot=3, s_tot=0, j_tot=3, f_tot=3.5, species="Yb171"), # "6snf 1F3" + AngularKetLS(l_c=0, l_r=3, l_tot=3, s_tot=1, j_tot=3, f_tot=3.5, species="Yb171"), # "6snf 3F3" + AngularKetDummy("4f13 5d 6snl a", f_tot=3.5), + AngularKetDummy("4f13 5d 6snl b", f_tot=3.5), + AngularKetDummy("4f13 5d 6snl c", f_tot=3.5), + AngularKetDummy("4f13 5d 6snl d", f_tot=3.5), + AngularKetDummy("4f13 5d 6snl e", f_tot=3.5), + AngularKetLS(l_c=0, l_r=3, l_tot=3, s_tot=1, j_tot=4, f_tot=3.5, species="Yb171"), # "6snf 3F4" + ] + outer_channels = [ + AngularKetFJ(l_c=0, l_r=3, j_c=0.5, f_c=1, j_r=3.5, f_tot=3.5, species="Yb171"), + AngularKetFJ(l_c=0, l_r=3, j_c=0.5, f_c=1, j_r=2.5, f_tot=3.5, species="Yb171"), + AngularKetDummy("4f13 5d 6snl a", f_tot=3.5), + AngularKetDummy("4f13 5d 6snl b", f_tot=3.5), + AngularKetDummy("4f13 5d 6snl c", f_tot=3.5), + AngularKetDummy("4f13 5d 6snl d", f_tot=3.5), + AngularKetDummy("4f13 5d 6snl e", f_tot=3.5), + AngularKetFJ(l_c=0, l_r=3, j_c=0.5, f_c=0, j_r=3.5, f_tot=3.5, species="Yb171"), + ] + + eigen_quantum_defects = [ + [0.277086649, -13.290196, 0], + [0.0719837014, -0.754736076, 0], + [0.251457795, 0, 0], + [0.227434828, 0, 0], + [0.175780645, 0, 0], + [0.196547521, 0, 0], + [0.21440857, 0, 0], + [0.0834193873, -1.11453386, -1545.71844], + ] + mixing_angles = [ + (0, 1, [0.0209955122, 0.251041249]), + (0, 2, -0.0585753224), + (0, 3, -0.0750574327), + (0, 4, 0.122671919), + (0, 5, -0.0401036164), + (0, 6, 0.0654271994), + (1, 2, -0.0683007974), + (1, 3, 0.035415976), + (1, 4, -0.0327625807), + (1, 5, -0.050225071), + (1, 6, 0.0455759316), + ] + + +class Yb171_F45_HighN(FModel): + species_name = "Yb171_mqdt" + name = "F F=4.5, nu > 20" + f_tot = 4.5 + nu_range = (20.0, np.inf) + reference = "10.1103/PhysRevA.112.042817" + + inner_channels = [ + AngularKetLS(l_c=0, l_r=3, l_tot=3, s_tot=1, j_tot=4, f_tot=4.5, species="Yb171"), # "6snf 3F4" + ] + outer_channels = [ + AngularKetFJ(l_c=0, l_r=3, j_c=0.5, f_c=1, j_r=3.5, f_tot=4.5, species="Yb171"), + ] + + eigen_quantum_defects = [ + [0.0834193873, -1.11453386, -1545.71844], + ] + mixing_angles = [] + + +class Yb171_G25_HighN(FModel): + species_name = "Yb171_mqdt" + name = "G F=2.5, nu > 25" + f_tot = 2.5 + nu_range = (25.0, np.inf) + reference = "10.1103/PhysRevA.112.042817" + + inner_channels = [ + AngularKetLS(l_c=0, l_r=4, l_tot=4, s_tot=1, j_tot=3, f_tot=2.5, species="Yb171"), # "6sng 3G3" + ] + outer_channels = [ + AngularKetFJ(l_c=0, l_r=4, j_c=0.5, f_c=1, j_r=3.5, f_tot=2.5, species="Yb171"), + ] + + eigen_quantum_defects = [ + [0.02613255, -0.14203905], + ] + mixing_angles = [] + + +class Yb171_G35_HighN(FModel): + species_name = "Yb171_mqdt" + name = "G F=3.5, nu > 25" + f_tot = 3.5 + nu_range = (25.0, np.inf) + reference = "10.1103/PhysRevA.112.042817" + + inner_channels = [ + AngularKetJJ(l_c=0, l_r=4, j_c=0.5, j_r=4.5, j_tot=4, f_tot=3.5, species="Yb171"), # "6sng +G4" + AngularKetJJ(l_c=0, l_r=4, j_c=0.5, j_r=3.5, j_tot=4, f_tot=3.5, species="Yb171"), # "6sng -G4" + AngularKetJJ(l_c=0, l_r=4, j_c=0.5, j_r=3.5, j_tot=3, f_tot=3.5, species="Yb171"), # "6sng 3G3" + ] + outer_channels = [ + AngularKetFJ(l_c=0, l_r=4, j_c=0.5, f_c=1, j_r=4.5, f_tot=3.5, species="Yb171"), + AngularKetFJ(l_c=0, l_r=4, j_c=0.5, f_c=1, j_r=3.5, f_tot=3.5, species="Yb171"), + AngularKetFJ(l_c=0, l_r=4, j_c=0.5, f_c=0, j_r=3.5, f_tot=3.5, species="Yb171"), + ] + + eigen_quantum_defects = [ + [0.02628545, -0.13182564], + [0.02548145, -0.12028462], + [0.02613255, -0.14203905], + ] + mixing_angles = [ + (0, 1, -0.089123698), + ] + + +class Yb171_G45_HighN(FModel): + species_name = "Yb171_mqdt" + name = "G F=4.5, nu > 25" + f_tot = 4.5 + nu_range = (25.0, np.inf) + reference = "10.1103/PhysRevA.112.042817" + + inner_channels = [ + AngularKetJJ(l_c=0, l_r=4, j_c=0.5, j_r=4.5, j_tot=4, f_tot=4.5, species="Yb171"), # "6sng +G4" + AngularKetJJ(l_c=0, l_r=4, j_c=0.5, j_r=3.5, j_tot=4, f_tot=4.5, species="Yb171"), # "6sng -G4" + AngularKetJJ(l_c=0, l_r=4, j_c=0.5, j_r=4.5, j_tot=5, f_tot=4.5, species="Yb171"), # "6sng 3G5" + ] + outer_channels = [ + AngularKetFJ(l_c=0, l_r=4, j_c=0.5, f_c=1, j_r=4.5, f_tot=4.5, species="Yb171"), + AngularKetFJ(l_c=0, l_r=4, j_c=0.5, f_c=1, j_r=3.5, f_tot=4.5, species="Yb171"), + AngularKetFJ(l_c=0, l_r=4, j_c=0.5, f_c=0, j_r=4.5, f_tot=4.5, species="Yb171"), + ] + + eigen_quantum_defects = [ + [0.02628545, -0.13182564], + [0.02548145, -0.12028462], + [0.02536571, -0.18507079], + ] + mixing_angles = [ + (0, 1, -0.089123698), + ] + + +class Yb171_G55_HighN(FModel): + species_name = "Yb171_mqdt" + name = "G F=5.5, nu > 25" + f_tot = 5.5 + nu_range = (25.0, np.inf) + reference = "10.1103/PhysRevA.112.042817" + + inner_channels = [ + AngularKetLS(l_c=0, l_r=4, l_tot=4, s_tot=1, j_tot=5, f_tot=5.5, species="Yb171"), # "6sng 3G5" + ] + outer_channels = [ + AngularKetFJ(l_c=0, l_r=4, j_c=0.5, f_c=1, j_r=4.5, f_tot=5.5, species="Yb171"), + ] + + eigen_quantum_defects = [ + [0.02536571, -0.18507079], + ] + mixing_angles = [] + + +# -------------------------------------------------------- +# MQDT models valid at small n +# -------------------------------------------------------- + + +class Yb171_S05_LowN(FModel): + species_name = "Yb171_mqdt" + name = "S F=0.5, 2 < nu < 26" + f_tot = 0.5 + nu_range = (2.0, 26.0) + reference = "10.1103/PhysRevX.15.011009" + + inner_channels = [ + AngularKetLS(l_c=0, l_r=0, l_tot=0, s_tot=0, j_tot=0, f_tot=0.5, species="Yb171"), # "6sns 1S0" + AngularKetDummy("4f13 5d 6snl a", f_tot=0.5), + AngularKetLS(l_c=1, l_r=1, l_tot=0, s_tot=0, j_tot=0, f_tot=0.5, species="Yb171"), # "6pnp 1S0" + AngularKetDummy("4f13 5d 6snl b", f_tot=0.5), + AngularKetLS(l_c=1, l_r=1, l_tot=1, s_tot=1, j_tot=0, f_tot=0.5, species="Yb171"), # "6pnp 3P0" + AngularKetDummy("4f13 5d 6snl c", f_tot=0.5), + AngularKetLS(l_c=0, l_r=0, l_tot=0, s_tot=1, j_tot=1, f_tot=0.5, species="Yb171"), # "6sns 3S1" + ] + outer_channels = [ + AngularKetFJ(l_c=0, l_r=0, j_c=0.5, f_c=0, j_r=0.5, f_tot=0.5, species="Yb171"), + AngularKetDummy("4f13 5d 6snl a", f_tot=0.5), + AngularKetFJ(l_c=1, l_r=1, j_c=1.5, f_c=1, j_r=1.5, f_tot=0.5, species="Yb171"), + AngularKetDummy("4f13 5d 6snl b", f_tot=0.5), + AngularKetFJ(l_c=1, l_r=1, j_c=0.5, f_c=Unknown, j_r=0.5, f_tot=0.5, species="Yb171"), + AngularKetDummy("4f13 5d 6snl c", f_tot=0.5), + AngularKetFJ(l_c=0, l_r=0, j_c=0.5, f_c=1, j_r=0.5, f_tot=0.5, species="Yb171"), + ] + + eigen_quantum_defects = [ + [0.357488757, 0.163255076, 0], + [0.203917828, 0, 0], + [0.116813499, 0, 0], + [0.287210377, 0, 0], + [0.247550262, 0, 0], + [0.148686263, 0, 0], + [0.432841, 0.724559, -1.95424], + ] + mixing_angles = [ + (0, 1, 0.13179534), + (0, 2, 0.29748039), + (0, 3, 0.0553920359), + (2, 3, 0.100843905), + (2, 4, 0.10317753), + (0, 5, 0.137709223), + ] + + manual_frame_transformation_outer_inner = np.array( + [ + [1 / 2, 0, 0, 0, 0, 0, np.sqrt(3) / 2], + [0, 1, 0, 0, 0, 0, 0], + [0, 0, np.sqrt(2 / 3), 0, -np.sqrt(1 / 3), 0, 0], + [0, 0, 0, 1, 0, 0, 0], + [0, 0, np.sqrt(1 / 3), 0, np.sqrt(2 / 3), 0, 0], + [0, 0, 0, 0, 0, 1, 0], + [np.sqrt(3) / 2, 0, 0, 0, 0, 0, -1 / 2], + ] + ) + + +class Yb171_S15_LowN(FModel): + species_name = "Yb171_mqdt" + name = "S F=1.5, 2 < nu < 26" + f_tot = 1.5 + nu_range = (2.0, 26.0) + reference = None + + inner_channels = [ + AngularKetLS(l_c=0, l_r=0, l_tot=0, s_tot=1, j_tot=1, f_tot=1.5, species="Yb171"), # "6sns 3S1" + ] + outer_channels = [ + AngularKetFJ(l_c=0, l_r=0, j_c=0.5, f_c=1, j_r=0.5, f_tot=1.5, species="Yb171"), + ] + + eigen_quantum_defects = [ + [0.432841, 0.724559, -1.95424], + ] + mixing_angles = [] + + +class Yb171_P05_Lowest(FModel): + species_name = "Yb171_mqdt" + name = "P F=0.5, 1.5 < nu < 2.5" + f_tot = 0.5 + nu_range = (1.5, 2.5) + reference = None + + inner_channels = [ + AngularKetLS(l_c=0, l_r=1, l_tot=1, s_tot=0, j_tot=1, f_tot=0.5, species="Yb171"), # "6snp 1P1" + AngularKetLS(l_c=0, l_r=1, l_tot=1, s_tot=1, j_tot=1, f_tot=0.5, species="Yb171"), # "6snp 3P1" + AngularKetLS(l_c=0, l_r=1, l_tot=1, s_tot=1, j_tot=0, f_tot=0.5, species="Yb171"), # "6snp 3P0" + ] + outer_channels = [ + AngularKetFJ(l_c=0, l_r=1, j_c=0.5, f_c=1, j_r=1.5, f_tot=0.5, species="Yb171"), + AngularKetFJ(l_c=0, l_r=1, j_c=0.5, f_c=1, j_r=0.5, f_tot=0.5, species="Yb171"), + AngularKetFJ(l_c=0, l_r=1, j_c=0.5, f_c=0, j_r=0.5, f_tot=0.5, species="Yb171"), + ] + + eigen_quantum_defects = [ + [0.161083, 0], + [0.920424, 0], + [0.180701, 0], + ] + mixing_angles = [ + (0, 1, [-0.426128, 6.272986]), + ] + + +class Yb171_P05_LowN(FModel): + species_name = "Yb171_mqdt" + name = "P F=0.5, 2.9 < nu < 5.9" + f_tot = 0.5 + nu_range = (2.9, 5.9) + reference = None + + inner_channels = [ + AngularKetLS(l_c=0, l_r=1, l_tot=1, s_tot=0, j_tot=1, f_tot=0.5, species="Yb171"), # "6snp 1P1" + AngularKetLS(l_c=0, l_r=1, l_tot=1, s_tot=1, j_tot=1, f_tot=0.5, species="Yb171"), # "6snp 3P1" + AngularKetLS(l_c=0, l_r=1, l_tot=1, s_tot=1, j_tot=0, f_tot=0.5, species="Yb171"), # "6snp 3P0" + ] + outer_channels = [ + AngularKetFJ(l_c=0, l_r=1, j_c=0.5, f_c=1, j_r=1.5, f_tot=0.5, species="Yb171"), + AngularKetFJ(l_c=0, l_r=1, j_c=0.5, f_c=1, j_r=0.5, f_tot=0.5, species="Yb171"), + AngularKetFJ(l_c=0, l_r=1, j_c=0.5, f_c=0, j_r=0.5, f_tot=0.5, species="Yb171"), + ] + + eigen_quantum_defects = [ + [0.967223, -3.03997, 0.569205], + [0.967918, 0.25116, 0.868505], + [0.969279, 0.288219, 1.36228], + ] + mixing_angles = [] + + +class Yb171_P15_Lowest(FModel): + species_name = "Yb171_mqdt" + name = "P F=1.5, 1.5 < nu < 2.5" + f_tot = 1.5 + nu_range = (1.5, 2.5) + reference = None + + inner_channels = [ + AngularKetLS(l_c=0, l_r=1, l_tot=1, s_tot=0, j_tot=1, f_tot=1.5, species="Yb171"), # "6snp 1P1" + AngularKetLS(l_c=0, l_r=1, l_tot=1, s_tot=1, j_tot=1, f_tot=1.5, species="Yb171"), # "6snp 3P1" + AngularKetLS(l_c=0, l_r=1, l_tot=1, s_tot=1, j_tot=2, f_tot=1.5, species="Yb171"), # "6snp 3P2" + ] + outer_channels = [ + AngularKetFJ(l_c=0, l_r=1, j_c=0.5, f_c=1, j_r=1.5, f_tot=1.5, species="Yb171"), + AngularKetFJ(l_c=0, l_r=1, j_c=0.5, f_c=1, j_r=0.5, f_tot=1.5, species="Yb171"), + AngularKetFJ(l_c=0, l_r=1, j_c=0.5, f_c=0, j_r=1.5, f_tot=1.5, species="Yb171"), + ] + + eigen_quantum_defects = [ + [0.161083, 0], + [0.920424, 0], + [0.110501, 0], + ] + mixing_angles = [ + (0, 1, [-0.426128, 6.272986]), + ] + + +class Yb171_P15_LowN(FModel): + species_name = "Yb171_mqdt" + name = "P F=1.5, 3 < nu < 10" + f_tot = 1.5 + nu_range = (3.0, 10.0) + reference = None + + inner_channels = [ + AngularKetLS(l_c=0, l_r=1, l_tot=1, s_tot=0, j_tot=1, f_tot=1.5, species="Yb171"), # "6snp 1P1" + AngularKetLS(l_c=0, l_r=1, l_tot=1, s_tot=1, j_tot=1, f_tot=1.5, species="Yb171"), # "6snp 3P1" + AngularKetDummy("4f13 5d 6snl a", f_tot=1.5), + AngularKetDummy("4f13 5d 6snl b", f_tot=1.5), + AngularKetDummy("4f13 5d 6snl c", f_tot=1.5), + AngularKetDummy("4f13 5d 6snl d", f_tot=1.5), + AngularKetLS(l_c=0, l_r=1, l_tot=1, s_tot=1, j_tot=2, f_tot=1.5, species="Yb171"), # "6snp 3P2" + AngularKetDummy("4f13 5d 6snl e", f_tot=1.5), + AngularKetDummy("4f13 5d 6snl f", f_tot=1.5), + AngularKetDummy("4f13 5d 6snl g", f_tot=1.5), + ] + outer_channels = [ + AngularKetFJ(l_c=0, l_r=1, j_c=0.5, f_c=1, j_r=1.5, f_tot=1.5, species="Yb171"), + AngularKetFJ(l_c=0, l_r=1, j_c=0.5, f_c=1, j_r=0.5, f_tot=1.5, species="Yb171"), + AngularKetDummy("4f13 5d 6snl a", f_tot=1.5), + AngularKetDummy("4f13 5d 6snl b", f_tot=1.5), + AngularKetDummy("4f13 5d 6snl c", f_tot=1.5), + AngularKetDummy("4f13 5d 6snl d", f_tot=1.5), + AngularKetFJ(l_c=0, l_r=1, j_c=0.5, f_c=0, j_r=1.5, f_tot=1.5, species="Yb171"), + AngularKetDummy("4f13 5d 6snl e", f_tot=1.5), + AngularKetDummy("4f13 5d 6snl f", f_tot=1.5), + AngularKetDummy("4f13 5d 6snl g", f_tot=1.5), + ] + + eigen_quantum_defects = [ + [0.967223, -3.03997, 0.569205], + [0.967918, 0.25116, 0.868505], + [0.228828720, 0, 0], + [0.205484818, 0, 0], + [0.193528629, 0, 0], + [0.181385000, 0, 0], + [0.906105, 0.383471, 1.23512], + [0.236866903, 0, 0], + [0.221055883, 0, 0], + [0.185599376, 0, 0], + ] + mixing_angles = [ + (0, 1, [-0.08410871, 120.37555, -9314.23]), + (0, 2, -0.073904060), + (0, 3, -0.063632668), + (0, 4, -0.021924569), + (0, 5, -0.106678810), + (1, 2, 0.032556999), + (1, 3, 0.054105142), + (1, 4, 0.086127672), + (1, 5, 0.053804487), + (6, 7, 0.071426685), + (6, 8, 0.027464110), + (6, 9, -0.029741862), + ] + + +class Yb171_P25_Lowest(FModel): + species_name = "Yb171_mqdt" + name = "P F=2.5, 1.5 < nu < 4.5" + f_tot = 2.5 + nu_range = (1.5, 4.5) + reference = None + + inner_channels = [ + AngularKetLS(l_c=0, l_r=1, l_tot=1, s_tot=1, j_tot=2, f_tot=2.5, species="Yb171"), # "6snp 3P2" + ] + outer_channels = [ + AngularKetFJ(l_c=0, l_r=1, j_c=0.5, f_c=1, j_r=1.5, f_tot=2.5, species="Yb171"), + ] + + eigen_quantum_defects = [ + [0.906105, 0.383471, 1.23512], + ] + mixing_angles = [] + + +class Yb171_P25_LowN(FModel): + species_name = "Yb171_mqdt" + name = "P F=2.5, 5 < nu < 20" + f_tot = 2.5 + nu_range = (5.0, 20.0) + reference = "10.1103/PhysRevX.15.011009" + + inner_channels = [ + AngularKetLS(l_c=0, l_r=1, l_tot=1, s_tot=1, j_tot=2, f_tot=2.5, species="Yb171"), # "6snp 3P2" + AngularKetDummy("4f13 5d 6snl a", f_tot=2.5), + AngularKetDummy("4f13 5d 6snl b", f_tot=2.5), + AngularKetDummy("4f13 5d 6snl c", f_tot=2.5), + ] + outer_channels = [ + AngularKetFJ(l_c=0, l_r=1, j_c=0.5, f_c=1, j_r=1.5, f_tot=2.5, species="Yb171"), + AngularKetDummy("4f13 5d 6snl a", f_tot=2.5), + AngularKetDummy("4f13 5d 6snl b", f_tot=2.5), + AngularKetDummy("4f13 5d 6snl c", f_tot=2.5), + ] + + eigen_quantum_defects = [ + [0.925121305, -2.73247165, 74.664989], + [0.230133261, 0, 0], + [0.209638118, 0, 0], + [0.186228192, 0, 0], + ] + mixing_angles = [ + (0, 1, 0.0706666127), + (0, 2, 0.0232711158), + (0, 3, -0.0292153659), + ] + + +class Yb171_D05_LowN(FModel): + species_name = "Yb171_mqdt" + name = "D F=0.5, 2 < nu < 30" + f_tot = 0.5 + nu_range = (2.0, 30.0) + reference = None + + inner_channels = [ + AngularKetLS(l_c=0, l_r=2, l_tot=2, s_tot=1, j_tot=1, f_tot=0.5, species="Yb171"), # "6snd 3D1" + ] + outer_channels = [ + AngularKetFJ(l_c=0, l_r=2, j_c=0.5, f_c=1, j_r=1.5, f_tot=0.5, species="Yb171"), + ] + + eigen_quantum_defects = [ + [0.758222, -0.017906, 3.392161], + ] + mixing_angles = [] + + +class Yb171_D15_LowN(FModel): + species_name = "Yb171_mqdt" + name = "D F=1.5, 2 < nu < 30" + f_tot = 1.5 + nu_range = (2.0, 30.0) + reference = "arXiv:2507.11487v1" + + inner_channels = [ + AngularKetLS(l_c=0, l_r=2, l_tot=2, s_tot=0, j_tot=2, f_tot=1.5, species="Yb171"), # "6snd 1D2" + AngularKetLS(l_c=0, l_r=2, l_tot=2, s_tot=1, j_tot=2, f_tot=1.5, species="Yb171"), # "6snd 3D2" + AngularKetDummy("4f13 5d 6snl a", f_tot=1.5), + AngularKetDummy("4f13 5d 6snl b", f_tot=1.5), + AngularKetLS(l_c=1, l_r=1, l_tot=2, s_tot=0, j_tot=2, f_tot=1.5, species="Yb171"), # "6pnp 1D2" + AngularKetLS(l_c=0, l_r=2, l_tot=2, s_tot=1, j_tot=1, f_tot=1.5, species="Yb171"), # "6snd 3D1" + ] + outer_channels = [ + AngularKetFJ(l_c=0, l_r=2, j_c=0.5, f_c=1, j_r=2.5, f_tot=1.5, species="Yb171"), + AngularKetFJ(l_c=0, l_r=2, j_c=0.5, f_c=1, j_r=1.5, f_tot=1.5, species="Yb171"), + AngularKetDummy("4f13 5d 6snl a", f_tot=1.5), + AngularKetDummy("4f13 5d 6snl b", f_tot=1.5), + AngularKetFJ(l_c=1, l_r=1, j_c=Unknown, f_c=Unknown, j_r=Unknown, f_tot=1.5, species="Yb171"), + AngularKetFJ(l_c=0, l_r=2, j_c=0.5, f_c=0, j_r=1.5, f_tot=1.5, species="Yb171"), + ] + manual_frame_transformation_outer_inner = np.array( + [ + [-np.sqrt(3 / 5), -np.sqrt(2 / 5), 0, 0, 0, 0], + [np.sqrt(3 / 5) / 2, -3 / (2 * np.sqrt(10)), 0, 0, 0, np.sqrt(5 / 2) / 2], + [0, 0, 1, 0, 0, 0], + [0, 0, 0, 1, 0, 0], + [0, 0, 0, 0, 1, 0], + [-1 / 2, np.sqrt(3 / 2) / 2, 0, 0, 0, np.sqrt(3 / 2) / 2], + ] + ) + + eigen_quantum_defects = [ + [0.730541589, -0.0967938662, 0], + [0.751542685, 0.00038836127, 0], + [0.195864083, 0, 0], + [0.235944408, 0, 0], + [0.147483609, 0, 0], + [0.758222, -0.017906, 3.392161], + ] + mixing_angles = [ + (0, 1, 0.220048245), + (0, 2, 0.00427599), + (0, 3, 0.0381563093), + (1, 3, -0.00700797918), + (0, 4, 0.109380331), + (1, 4, 0.0635544456), + ] + + +class Yb171_D25_LowN(FModel): + species_name = "Yb171_mqdt" + name = "D F=2.5, 2 < nu < 30" + f_tot = 2.5 + nu_range = (2.0, 30.0) + reference = "arXiv:2507.11487v1" + + inner_channels = [ + AngularKetLS(l_c=0, l_r=2, l_tot=2, s_tot=0, j_tot=2, f_tot=2.5, species="Yb171"), # "6snd 1D2" + AngularKetLS(l_c=0, l_r=2, l_tot=2, s_tot=1, j_tot=2, f_tot=2.5, species="Yb171"), # "6snd 3D2" + AngularKetDummy("4f13 5d 6snl a", f_tot=2.5), + AngularKetDummy("4f13 5d 6snl b", f_tot=2.5), + AngularKetLS(l_c=1, l_r=1, l_tot=2, s_tot=0, j_tot=2, f_tot=2.5, species="Yb171"), # "6pnp 1D2" + AngularKetLS(l_c=0, l_r=2, l_tot=2, s_tot=1, j_tot=3, f_tot=2.5, species="Yb171"), # "6snd 3D3" + ] + outer_channels = [ + AngularKetFJ(l_c=0, l_r=2, j_c=0.5, f_c=1, j_r=2.5, f_tot=2.5, species="Yb171"), + AngularKetFJ(l_c=0, l_r=2, j_c=0.5, f_c=1, j_r=1.5, f_tot=2.5, species="Yb171"), + AngularKetDummy("4f13 5d 6snl a", f_tot=2.5), + AngularKetDummy("4f13 5d 6snl b", f_tot=2.5), + AngularKetFJ(l_c=1, l_r=1, j_c=Unknown, f_c=Unknown, j_r=Unknown, f_tot=2.5, species="Yb171"), + AngularKetFJ(l_c=0, l_r=2, j_c=0.5, f_c=0, j_r=2.5, f_tot=2.5, species="Yb171"), + ] + manual_frame_transformation_outer_inner = np.array( + [ + [np.sqrt(7 / 5) / 2, np.sqrt(7 / 30), 0, 0, 0, -np.sqrt(5 / 3) / 2], + [-np.sqrt(2 / 5), np.sqrt(3 / 5), 0, 0, 0, 0], + [0, 0, 1, 0, 0, 0], + [0, 0, 0, 1, 0, 0], + [0, 0, 0, 0, 1, 0], + [1 / 2, np.sqrt(1 / 6), 0, 0, 0, np.sqrt(7 / 3) / 2], + ] + ) + + eigen_quantum_defects = [ + [0.730541589, -0.0967938662, 0], + [0.751542685, 0.00038836127, 0], + [0.195864083, 0, 0], + [0.235944408, 0, 0], + [0.147483609, 0, 0], + [0.734512, -0.019501, 3.459114], + ] + mixing_angles = [ + (0, 1, [0.220048245, -14.9486]), + (0, 2, 0.00427599), + (0, 3, 0.0381563093), + (1, 3, -0.00700797918), + (0, 4, 0.109380331), + (1, 4, 0.0635544456), + ] + + +class Yb171_D35_LowN(FModel): + species_name = "Yb171_mqdt" + name = "D F=3.5, 2 < nu < 14" + f_tot = 3.5 + nu_range = (2.0, 14.0) + reference = None + + inner_channels = [ + AngularKetLS(l_c=0, l_r=2, l_tot=2, s_tot=1, j_tot=3, f_tot=3.5, species="Yb171"), # "6snd 3D3" + ] + outer_channels = [ + AngularKetFJ(l_c=0, l_r=2, j_c=0.5, f_c=1, j_r=2.5, f_tot=3.5, species="Yb171"), + ] + + eigen_quantum_defects = [ + [0.734512, -0.019501, 3.459114], + ] + mixing_angles = [] diff --git a/src/rydstate/species/mqdt/yb173.py b/src/rydstate/species/mqdt/yb173.py new file mode 100644 index 0000000..9a9af1e --- /dev/null +++ b/src/rydstate/species/mqdt/yb173.py @@ -0,0 +1,394 @@ +from __future__ import annotations + +import numpy as np + +from rydstate.angular.angular_ket import AngularKetFJ, AngularKetLS +from rydstate.angular.angular_ket_dummy import AngularKetDummy +from rydstate.angular.core_ket_base import CoreKet, CoreKetDummy +from rydstate.angular.utils import Unknown +from rydstate.species.mqdt.fmodel import FModel +from rydstate.species.mqdt.species_object_mqdt import SpeciesObjectMQDT +from rydstate.units import electron_mass, rydberg_constant + + +class Ytterbium173MQDT(SpeciesObjectMQDT): + name = "Yb173_mqdt" + Z = 70 + i_c = 2.5 + number_valence_electrons = 2 + + potential_type_default = "model_potential_fei_2009" + + # https://iopscience.iop.org/article/10.1088/1674-1056/18/10/025 + model_potential_parameter_fei_2009 = (0.8704, 22.0040, 0.1513, 0.3306) + + _isotope_mass = 172.938216212 # u + _corrected_rydberg_constant = ( + rydberg_constant.m / (1 + electron_mass.to("u").m / _isotope_mass), + None, + str(rydberg_constant.u), + ) + + _ionization_threshold_dict = { + CoreKet(i_c, 0.5, 0, 0.5, 2): (50443.291203, None, "1/cm"), + CoreKet(i_c, 0.5, 0, 0.5, 3): (50442.941262, None, "1/cm"), + CoreKet(i_c, 0.5, 1, 0.5, Unknown): (77504.98, None, "1/cm"), + CoreKet(i_c, 0.5, 1, 1.5, Unknown): (80835.39, None, "1/cm"), + CoreKetDummy("4f13 5d 6s"): (83967.7, None, "1/cm"), + } + core_ground_state = CoreKet(i_c, 0.5, 0, 0.5, 2) + nuclear_dipole = -0.68 + + +# -------------------------------------------------------- +# MQDT models valid at large n +# -------------------------------------------------------- + + +class Yb173_S15_HighN(FModel): + species_name = "Yb173_mqdt" + name = "S F=3/2, nu > 26" + f_tot = 1.5 + nu_range = (26.0, np.inf) + reference = "10.1103/PhysRevA.110.042821" + + inner_channels = [ + AngularKetLS(l_c=0, l_r=0, l_tot=0, s_tot=1, j_tot=1, f_tot=1.5, species="Yb173"), # "6sns 3S1" + ] + outer_channels = [ + AngularKetFJ(l_c=0, l_r=0, j_c=0.5, f_c=2, j_r=0.5, f_tot=1.5, species="Yb173"), + ] + + eigen_quantum_defects = [ + [0.438426851, 3.91762642, -10612.6828, 8017432.38, -2582622910.0], + ] + mixing_angles = [] + + +class Yb173_S25_HighN(FModel): + species_name = "Yb173_mqdt" + name = "S F=5/2, nu > 26" + f_tot = 2.5 + nu_range = (26.0, np.inf) + reference = "10.1103/PhysRevA.110.042821" + + inner_channels = [ + AngularKetLS(l_c=0, l_r=0, l_tot=0, s_tot=0, j_tot=0, f_tot=2.5, species="Yb173"), # "6sns 1S0" + AngularKetDummy("4f13 5d 6snl a", f_tot=2.5), + AngularKetLS(l_c=1, l_r=1, l_tot=0, s_tot=0, j_tot=0, f_tot=2.5, species="Yb173"), # "6pnp 1S0" + AngularKetDummy("4f13 5d 6snl b", f_tot=2.5), + AngularKetLS(l_c=1, l_r=1, l_tot=1, s_tot=1, j_tot=0, f_tot=2.5, species="Yb173"), # "6pnp 3P0" + AngularKetDummy("4f13 5d 6snl c", f_tot=2.5), + AngularKetLS(l_c=0, l_r=0, l_tot=0, s_tot=1, j_tot=1, f_tot=2.5, species="Yb173"), # "6sns 3S1" + ] + outer_channels = [ + AngularKetFJ(l_c=0, l_r=0, j_c=0.5, f_c=2, j_r=0.5, f_tot=2.5, species="Yb173"), + AngularKetDummy("4f13 5d 6snl a", f_tot=2.5), + AngularKetFJ(l_c=1, l_r=1, j_c=1.5, f_c=3, j_r=1.5, f_tot=2.5, species="Yb173"), + AngularKetDummy("4f13 5d 6snl b", f_tot=2.5), + AngularKetFJ(l_c=1, l_r=1, j_c=0.5, f_c=Unknown, j_r=0.5, f_tot=2.5, species="Yb173"), # Fc could be 2 or 3 + AngularKetDummy("4f13 5d 6snl c", f_tot=2.5), + AngularKetFJ(l_c=0, l_r=0, j_c=0.5, f_c=3, j_r=0.5, f_tot=2.5, species="Yb173"), + ] + manual_frame_transformation_outer_inner = np.array( + [ + [np.sqrt(5) / 2 / np.sqrt(3), 0, 0, 0, 0, 0, np.sqrt(7) / 2 / np.sqrt(3)], + [0, 1, 0, 0, 0, 0, 0], + [0, 0, -np.sqrt(2 / 3), 0, np.sqrt(1 / 3), 0, 0], + [0, 0, 0, 1, 0, 0, 0], + [0, 0, np.sqrt(1 / 3), 0, np.sqrt(2 / 3), 0, 0], + [0, 0, 0, 0, 0, 1, 0], + [np.sqrt(7) / 2 / np.sqrt(3), 0, 0, 0, 0, 0, -np.sqrt(5) / 2 / np.sqrt(3)], + ] + ) + + eigen_quantum_defects = [ + [0.357519763, 0.298712849, 0, 0, 0], + [0.203907536, 0, 0, 0, 0], + [0.116803536, 0, 0, 0, 0], + [0.286731074, 0, 0, 0, 0], + [0.248113946, 0, 0, 0, 0], + [0.148678953, 0, 0, 0, 0], + [0.438426851, 3.91762642, -10612.6828, 8017432.38, -2582622910.0], + ] + mixing_angles = [ + (0, 1, 0.131810463), + (0, 2, 0.297612147), + (0, 3, 0.055508821), + (2, 3, 0.101030515), + (2, 4, 0.102911159), + (0, 5, 0.137723736), + ] + + +class Yb173_S35_HighN(FModel): + species_name = "Yb173_mqdt" + name = "S F=7/2, nu > 26" + f_tot = 3.5 + nu_range = (26.0, np.inf) + reference = "10.1103/PhysRevA.110.042821" + + inner_channels = [ + AngularKetLS(l_c=0, l_r=0, l_tot=0, s_tot=1, j_tot=1, f_tot=3.5, species="Yb173"), # "6sns 3S1" + ] + outer_channels = [ + AngularKetFJ(l_c=0, l_r=0, j_c=0.5, f_c=3, j_r=0.5, f_tot=3.5, species="Yb173"), + ] + + eigen_quantum_defects = [ + [0.438426851, 3.91762642, -10612.6828, 8017432.38, -2582622910.0], + ] + mixing_angles = [] + + +class Yb173_P05_HighN(FModel): + species_name = "Yb173_mqdt" + name = "P F=1/2, nu > 10" + f_tot = 0.5 + nu_range = (10.0, np.inf) + reference = "10.1103/PhysRevA.110.042821" + + inner_channels = [ + AngularKetLS(l_c=0, l_r=1, l_tot=1, s_tot=1, j_tot=2, f_tot=0.5, species="Yb173"), # "6snp 3P2" + AngularKetDummy("4f13 5d 6snl a", f_tot=0.5), + AngularKetDummy("4f13 5d 6snl b", f_tot=0.5), + AngularKetDummy("4f13 5d 6snl c", f_tot=0.5), + ] + outer_channels = [ + AngularKetFJ(l_c=0, l_r=1, j_c=0.5, f_c=2, j_r=1.5, f_tot=0.5, species="Yb173"), + AngularKetDummy("4f13 5d 6snl a", f_tot=0.5), + AngularKetDummy("4f13 5d 6snl b", f_tot=0.5), + AngularKetDummy("4f13 5d 6snl c", f_tot=0.5), + ] + + eigen_quantum_defects = [ + [0.924825736, -3.542481644, 81.5334687], + [0.236866903, 0, 0], + [0.221055883, 0, 0], + [0.185599376, 0, 0], + ] + mixing_angles = [ + (0, 1, 0.071426685), + (0, 2, 0.027464110), + (0, 3, -0.029741862), + ] + + +class Yb173_P15_HighN(FModel): + species_name = "Yb173_mqdt" + name = "P F=3/2, nu > 10" + f_tot = 1.5 + nu_range = (10.0, np.inf) + reference = "10.1103/PhysRevA.110.042821" + + inner_channels = [ + AngularKetLS(l_c=0, l_r=1, l_tot=1, s_tot=0, j_tot=1, f_tot=1.5, species="Yb173"), # "6snp 1P1" + AngularKetLS(l_c=0, l_r=1, l_tot=1, s_tot=1, j_tot=1, f_tot=1.5, species="Yb173"), # "6snp 3P1" + AngularKetDummy("4f13 5d 6snl a", f_tot=1.5), + AngularKetDummy("4f13 5d 6snl b", f_tot=1.5), + AngularKetDummy("4f13 5d 6snl c", f_tot=1.5), + AngularKetDummy("4f13 5d 6snl d", f_tot=1.5), + AngularKetLS(l_c=0, l_r=1, l_tot=1, s_tot=1, j_tot=2, f_tot=1.5, species="Yb173"), # "6snp 3P2" + AngularKetDummy("4f13 5d 6snl e", f_tot=1.5), + AngularKetDummy("4f13 5d 6snl f", f_tot=1.5), + AngularKetDummy("4f13 5d 6snl g", f_tot=1.5), + ] + outer_channels = [ + AngularKetFJ(l_c=0, l_r=1, j_c=0.5, f_c=2, j_r=0.5, f_tot=1.5, species="Yb173"), + AngularKetFJ(l_c=0, l_r=1, j_c=0.5, f_c=2, j_r=1.5, f_tot=1.5, species="Yb173"), + AngularKetDummy("4f13 5d 6snl a", f_tot=1.5), + AngularKetDummy("4f13 5d 6snl b", f_tot=1.5), + AngularKetDummy("4f13 5d 6snl c", f_tot=1.5), + AngularKetDummy("4f13 5d 6snl d", f_tot=1.5), + AngularKetFJ(l_c=0, l_r=1, j_c=0.5, f_c=3, j_r=1.5, f_tot=1.5, species="Yb173"), + AngularKetDummy("4f13 5d 6snl e", f_tot=1.5), + AngularKetDummy("4f13 5d 6snl f", f_tot=1.5), + AngularKetDummy("4f13 5d 6snl g", f_tot=1.5), + ] + + eigen_quantum_defects = [ + [0.921706585, 2.56569459, 0], + [0.979638580, -5.239904224, 0], + [0.228828720, 0, 0], + [0.205484818, 0, 0], + [0.193528629, 0, 0], + [0.181385000, 0, 0], + [0.924825736, -3.542481644, 81.5334687], + [0.236866903, 0, 0], + [0.221055883, 0, 0], + [0.185599376, 0, 0], + ] + mixing_angles = [ + (0, 1, [-0.087127227, 135.400009, -12985.0162]), + (0, 2, -0.073904060), + (0, 3, -0.063632668), + (0, 4, -0.021924569), + (0, 5, -0.106678810), + (1, 2, 0.032556999), + (1, 3, 0.054105142), + (1, 4, 0.086127672), + (1, 5, 0.053804487), + (6, 7, 0.071426685), + (6, 8, 0.027464110), + (6, 9, -0.029741862), + ] + + +class Yb173_P25_HighN(FModel): + species_name = "Yb173_mqdt" + name = "P F=5/2, nu > 10" + f_tot = 2.5 + nu_range = (10.0, np.inf) + reference = "10.1103/PhysRevA.110.042821" + + inner_channels = [ + AngularKetLS(l_c=0, l_r=1, l_tot=1, s_tot=0, j_tot=1, f_tot=2.5, species="Yb173"), # "6snp 1P1" + AngularKetLS(l_c=0, l_r=1, l_tot=1, s_tot=1, j_tot=1, f_tot=2.5, species="Yb173"), # "6snp 3P1" + AngularKetDummy("4f13 5d 6snl a", f_tot=2.5), + AngularKetDummy("4f13 5d 6snl b", f_tot=2.5), + AngularKetDummy("4f13 5d 6snl c", f_tot=2.5), + AngularKetDummy("4f13 5d 6snl d", f_tot=2.5), + AngularKetLS(l_c=0, l_r=1, l_tot=1, s_tot=1, j_tot=0, f_tot=2.5, species="Yb173"), # "6snp 3P0" + AngularKetDummy("4f13 5d 6snl e", f_tot=2.5), + AngularKetLS(l_c=0, l_r=1, l_tot=1, s_tot=1, j_tot=2, f_tot=2.5, species="Yb173"), # "6snp 3P2" + AngularKetDummy("4f13 5d 6snl f", f_tot=2.5), + AngularKetDummy("4f13 5d 6snl g", f_tot=2.5), + AngularKetDummy("4f13 5d 6snl h", f_tot=2.5), + ] + outer_channels = [ + AngularKetFJ(l_c=0, l_r=1, j_c=0.5, f_c=2, j_r=0.5, f_tot=2.5, species="Yb173"), + AngularKetFJ(l_c=0, l_r=1, j_c=0.5, f_c=2, j_r=1.5, f_tot=2.5, species="Yb173"), + AngularKetDummy("4f13 5d 6snl a", f_tot=2.5), + AngularKetDummy("4f13 5d 6snl b", f_tot=2.5), + AngularKetDummy("4f13 5d 6snl c", f_tot=2.5), + AngularKetDummy("4f13 5d 6snl d", f_tot=2.5), + AngularKetFJ(l_c=0, l_r=1, j_c=0.5, f_c=3, j_r=0.5, f_tot=2.5, species="Yb173"), + AngularKetDummy("4f13 5d 6snl e", f_tot=2.5), + AngularKetFJ(l_c=0, l_r=1, j_c=0.5, f_c=3, j_r=1.5, f_tot=2.5, species="Yb173"), + AngularKetDummy("4f13 5d 6snl f", f_tot=2.5), + AngularKetDummy("4f13 5d 6snl g", f_tot=2.5), + AngularKetDummy("4f13 5d 6snl h", f_tot=2.5), + ] + + eigen_quantum_defects = [ + [0.921706585, 2.56569459, 0], + [0.979638580, -5.239904224, 0], + [0.228828720, 0, 0], + [0.205484818, 0, 0], + [0.193528629, 0, 0], + [0.181385000, 0, 0], + [0.95356884, -0.28602498, 0], + [0.19845903, 0, 0], + [0.924825736, -3.542481644, 81.5334687], + [0.236866903, 0, 0], + [0.221055883, 0, 0], + [0.185599376, 0, 0], + ] + mixing_angles = [ + (0, 1, [-0.087127227, 135.400009, -12985.0162]), + (0, 2, -0.073904060), + (0, 3, -0.063632668), + (0, 4, -0.021924569), + (0, 5, -0.106678810), + (1, 2, 0.032556999), + (1, 3, 0.054105142), + (1, 4, 0.086127672), + (1, 5, 0.053804487), + (6, 7, 0.16328854), + (9, 8, -0.071426685), + (10, 8, -0.027464110), + (11, 8, 0.029741862), + ] + + +class Yb173_P35_HighN(FModel): + species_name = "Yb173_mqdt" + name = "P F=7/2, nu > 10" + f_tot = 3.5 + nu_range = (10.0, np.inf) + reference = "10.1103/PhysRevA.110.042821" + + inner_channels = [ + AngularKetLS(l_c=0, l_r=1, l_tot=1, s_tot=0, j_tot=1, f_tot=3.5, species="Yb173"), # "6snp 1P1" + AngularKetLS(l_c=0, l_r=1, l_tot=1, s_tot=1, j_tot=1, f_tot=3.5, species="Yb173"), # "6snp 3P1" + AngularKetDummy("4f13 5d 6snl a", f_tot=3.5), + AngularKetDummy("4f13 5d 6snl b", f_tot=3.5), + AngularKetDummy("4f13 5d 6snl c", f_tot=3.5), + AngularKetDummy("4f13 5d 6snl d", f_tot=3.5), + AngularKetLS(l_c=0, l_r=1, l_tot=1, s_tot=1, j_tot=2, f_tot=3.5, species="Yb173"), # "6snp 3P2" + AngularKetDummy("4f13 5d 6snl e", f_tot=3.5), + AngularKetDummy("4f13 5d 6snl f", f_tot=3.5), + AngularKetDummy("4f13 5d 6snl g", f_tot=3.5), + ] + outer_channels = [ + AngularKetFJ(l_c=0, l_r=1, j_c=0.5, f_c=2, j_r=1.5, f_tot=3.5, species="Yb173"), + AngularKetFJ(l_c=0, l_r=1, j_c=0.5, f_c=3, j_r=0.5, f_tot=3.5, species="Yb173"), + AngularKetDummy("4f13 5d 6snl a", f_tot=3.5), + AngularKetDummy("4f13 5d 6snl b", f_tot=3.5), + AngularKetDummy("4f13 5d 6snl c", f_tot=3.5), + AngularKetDummy("4f13 5d 6snl d", f_tot=3.5), + AngularKetFJ(l_c=0, l_r=1, j_c=0.5, f_c=3, j_r=1.5, f_tot=3.5, species="Yb173"), + AngularKetDummy("4f13 5d 6snl e", f_tot=3.5), + AngularKetDummy("4f13 5d 6snl f", f_tot=3.5), + AngularKetDummy("4f13 5d 6snl g", f_tot=3.5), + ] + + eigen_quantum_defects = [ + [0.921706585, 2.56569459, 0], + [0.979638580, -5.239904224, 0], + [0.228828720, 0, 0], + [0.205484818, 0, 0], + [0.193528629, 0, 0], + [0.181385000, 0, 0], + [0.924825736, -3.542481644, 81.5334687], + [0.236866903, 0, 0], + [0.221055883, 0, 0], + [0.185599376, 0, 0], + ] + mixing_angles = [ + (0, 1, [-0.087127227, 135.400009, -12985.0162]), + (0, 2, -0.073904060), + (0, 3, -0.063632668), + (0, 4, -0.021924569), + (0, 5, -0.106678810), + (1, 2, 0.032556999), + (1, 3, 0.054105142), + (1, 4, 0.086127672), + (1, 5, 0.053804487), + (6, 7, 0.071426685), + (6, 8, 0.027464110), + (6, 9, -0.029741862), + ] + + +class Yb173_P45_HighN(FModel): + species_name = "Yb173_mqdt" + name = "P F=9/2, nu > 10" + f_tot = 4.5 + nu_range = (10.0, np.inf) + reference = "10.1103/PhysRevA.110.042821" + + inner_channels = [ + AngularKetLS(l_c=0, l_r=1, l_tot=1, s_tot=1, j_tot=2, f_tot=4.5, species="Yb173"), # "6snp 3P2" + AngularKetDummy("4f13 5d 6snl a", f_tot=4.5), + AngularKetDummy("4f13 5d 6snl b", f_tot=4.5), + AngularKetDummy("4f13 5d 6snl c", f_tot=4.5), + ] + outer_channels = [ + AngularKetFJ(l_c=0, l_r=1, j_c=0.5, f_c=3, j_r=1.5, f_tot=4.5, species="Yb173"), + AngularKetDummy("4f13 5d 6snl a", f_tot=4.5), + AngularKetDummy("4f13 5d 6snl b", f_tot=4.5), + AngularKetDummy("4f13 5d 6snl c", f_tot=4.5), + ] + + eigen_quantum_defects = [ + [0.924825736, -3.542481644, 81.5334687], + [0.236866903, 0, 0], + [0.221055883, 0, 0], + [0.185599376, 0, 0], + ] + mixing_angles = [ + (0, 1, 0.071426685), + (0, 2, 0.027464110), + (0, 3, -0.029741862), + ] diff --git a/src/rydstate/species/mqdt/yb174.py b/src/rydstate/species/mqdt/yb174.py new file mode 100644 index 0000000..adcd416 --- /dev/null +++ b/src/rydstate/species/mqdt/yb174.py @@ -0,0 +1,652 @@ +from __future__ import annotations + +import numpy as np + +from rydstate.angular.angular_ket import AngularKetFJ, AngularKetJJ, AngularKetLS +from rydstate.angular.angular_ket_dummy import AngularKetDummy +from rydstate.angular.core_ket_base import CoreKet, CoreKetDummy +from rydstate.angular.utils import Unknown +from rydstate.species.mqdt.fmodel import FModel +from rydstate.species.mqdt.species_object_mqdt import SpeciesObjectMQDT +from rydstate.units import electron_mass, rydberg_constant + + +class Ytterbium174MQDT(SpeciesObjectMQDT): + name = "Yb174_mqdt" + Z = 70 + i_c = 0 + number_valence_electrons = 2 + + potential_type_default = "model_potential_fei_2009" + + # https://iopscience.iop.org/article/10.1088/1674-1056/18/10/025 + model_potential_parameter_fei_2009 = (0.8704, 22.0040, 0.1513, 0.3306) + + # https://physics.nist.gov/PhysRefData/Handbook/Tables/ytterbiumtable1.htm + _isotope_mass = 173.938859 # u + _corrected_rydberg_constant = ( + rydberg_constant.m / (1 + electron_mass.to("u").m / _isotope_mass), + None, + str(rydberg_constant.u), + ) + + _ionization_threshold_dict = { + CoreKet(i_c, 0.5, 0, 0.5): (50443.070393, None, "1/cm"), + CoreKet(i_c, 0.5, 1, 0.5): (77504.98, None, "1/cm"), + CoreKet(i_c, 0.5, 1, 1.5): (80835.39, None, "1/cm"), + CoreKet(i_c, 0.5, 1, Unknown): (79725.35, None, "1/cm"), + CoreKetDummy("4f13 5d 6s"): (83967.7, None, "1/cm"), + } + core_ground_state = CoreKet(i_c, 0.5, 0, 0.5) + nuclear_dipole = 2.1 + + +# -------------------------------------------------------- +# MQDT models valid at large n +# -------------------------------------------------------- + + +class Yb174_S0_HighN(FModel): + species_name = "Yb174_mqdt" + name = "S J=0, nu > 2" + f_tot = 0 + nu_range = (2.0, np.inf) + reference = "10.1103/PhysRevX.15.011009" + + inner_channels = [ + AngularKetLS(l_c=0, l_r=0, l_tot=0, s_tot=0, j_tot=0, species="Yb174"), # "6sns 1S0" + AngularKetDummy("4f13 5d 6snl a", f_tot=0), + AngularKetLS(l_c=1, l_r=1, l_tot=0, s_tot=0, j_tot=0, species="Yb174"), # "6pnp 1S0" + AngularKetDummy("4f13 5d 6snl b", f_tot=0), + AngularKetLS(l_c=1, l_r=1, l_tot=1, s_tot=1, j_tot=0, species="Yb174"), # "6pnp 3P0" + AngularKetDummy("4f13 5d 6snl c", f_tot=0), + ] + outer_channels = [ + AngularKetFJ(l_c=0, l_r=0, j_c=0.5, j_r=0.5, f_tot=0, species="Yb174"), + AngularKetDummy("4f13 5d 6snl a", f_tot=0), + AngularKetFJ(l_c=1, l_r=1, j_c=1.5, j_r=1.5, f_tot=0, species="Yb174"), + AngularKetDummy("4f13 5d 6snl b", f_tot=0), + AngularKetFJ(l_c=1, l_r=1, j_c=0.5, j_r=0.5, f_tot=0, species="Yb174"), + AngularKetDummy("4f13 5d 6snl c", f_tot=0), + ] + + eigen_quantum_defects = [ + [0.355101645, 0.277673956], + [0.204537535, 0], + [0.116393648, 0], + [0.295439966, 0], + [0.257664798, 0], + [0.155797119, 0], + ] + mixing_angles = [ + (0, 1, 0.126557575), + (0, 2, 0.300103593), + (0, 3, 0.056987912), + (2, 3, 0.114312578), + (2, 4, 0.0986363362), + (0, 5, 0.142498543), + ] + + +class Yb174_S1_HighN(FModel): + species_name = "Yb174_mqdt" + name = "S J=1, nu > 26" + f_tot = 1 + nu_range = (26.0, np.inf) + reference = "10.1103/PhysRevLett.128.033201" + + inner_channels = [ + AngularKetLS(l_c=0, l_r=0, l_tot=0, s_tot=1, j_tot=1, species="Yb174"), # "6sns 3S1" + ] + outer_channels = [ + AngularKetFJ(l_c=0, l_r=0, j_c=0.5, j_r=0.5, f_tot=1, species="Yb174"), + ] + + eigen_quantum_defects = [ + [0.4382, 4, -1e4, 8e6, -3e9], + ] + mixing_angles = [] + + +class Yb174_P0_HighN(FModel): + species_name = "Yb174_mqdt" + name = "P J=0, nu > 6" + f_tot = 0 + nu_range = (6.0, np.inf) + reference = "10.1103/PhysRevX.15.011009" + + inner_channels = [ + AngularKetLS(l_c=0, l_r=1, l_tot=1, s_tot=1, j_tot=0, species="Yb174"), # "6snp 3P0" + AngularKetDummy("4f13 5d 6snd", f_tot=0), + ] + outer_channels = [ + AngularKetFJ(l_c=0, l_r=1, j_c=0.5, j_r=0.5, f_tot=0, species="Yb174"), + AngularKetDummy("4f13 5d 6snd", f_tot=0), + ] + + eigen_quantum_defects = [ + [0.953661478, -0.287531374], + [0.198460766, 0], + ] + mixing_angles = [ + (0, 1, 0.163343232), + ] + + +class Yb174_P1_HighN(FModel): + species_name = "Yb174_mqdt" + name = "P J=1, nu > 6" + f_tot = 1 + nu_range = (6.0, np.inf) + reference = "10.1103/PhysRevX.15.011009" + + inner_channels = [ + AngularKetLS(l_c=0, l_r=1, l_tot=1, s_tot=0, j_tot=1, species="Yb174"), # "6snp 1P1" + AngularKetLS(l_c=0, l_r=1, l_tot=1, s_tot=1, j_tot=1, species="Yb174"), # "6snp 3P1" + AngularKetDummy("4f13 5d 6snl a", f_tot=1), + AngularKetDummy("4f13 5d 6snl b", f_tot=1), + AngularKetDummy("4f13 5d 6snl c", f_tot=1), + AngularKetDummy("4f13 5d 6snl d", f_tot=1), + ] + outer_channels = [ + AngularKetFJ(l_c=0, l_r=1, j_c=0.5, j_r=1.5, f_tot=1, species="Yb174"), + AngularKetFJ(l_c=0, l_r=1, j_c=0.5, j_r=0.5, f_tot=1, species="Yb174"), + AngularKetDummy("4f13 5d 6snl a", f_tot=1), + AngularKetDummy("4f13 5d 6snl b", f_tot=1), + AngularKetDummy("4f13 5d 6snl c", f_tot=1), + AngularKetDummy("4f13 5d 6snl d", f_tot=1), + ] + + eigen_quantum_defects = [ + [0.92271098, 2.6036257], + [0.98208719, -5.4562725], + [0.22851720, 0], + [0.20607759, 0], + [0.19352751, 0], + [0.18153094, 0], + ] + mixing_angles = [ + (0, 1, [-0.08410871, 120.37555, -9314.23]), + (0, 2, -0.07317986), + (0, 3, -0.06651879), + (0, 4, -0.02212194), + (0, 5, -0.10452109), + (1, 2, 0.02477464), + (1, 3, 0.05763934), + (1, 4, 0.0860644), + (1, 5, 0.04993818), + ] + + +class Yb174_P2_HighN(FModel): + species_name = "Yb174_mqdt" + name = "P J=2, nu > 5" + f_tot = 2 + nu_range = (5.0, np.inf) + reference = "10.1103/PhysRevX.15.011009" + + inner_channels = [ + AngularKetLS(l_c=0, l_r=1, l_tot=1, s_tot=1, j_tot=2, species="Yb174"), # "6snp 3P2" + AngularKetDummy("4f13 5d 6snl a", f_tot=2), + AngularKetDummy("4f13 5d 6snl b", f_tot=2), + AngularKetDummy("4f13 5d 6snl c", f_tot=2), + ] + outer_channels = [ + AngularKetFJ(l_c=0, l_r=1, j_c=0.5, j_r=1.5, f_tot=2, species="Yb174"), + AngularKetDummy("4f13 5d 6snl a", f_tot=2), + AngularKetDummy("4f13 5d 6snl b", f_tot=2), + AngularKetDummy("4f13 5d 6snl c", f_tot=2), + ] + + eigen_quantum_defects = [ + [0.925150932, -2.69197178, 66.7159709], + [0.230028034, 0, 0], + [0.209224174, 0, 0], + [0.186236574, 0, 0], + ] + mixing_angles = [ + (0, 1, 0.0706189664), + (0, 2, 0.0231221428), + (0, 3, -0.0291730345), + ] + + +class Yb174_D1_HighN(FModel): + species_name = "Yb174_mqdt" + name = "D J=1, nu > 26" + f_tot = 1 + nu_range = (26.0, np.inf) + reference = "10.1103/PhysRevX.15.011009" + + inner_channels = [ + AngularKetLS(l_c=0, l_r=2, l_tot=2, s_tot=1, j_tot=1, species="Yb174"), # "6snd 3D1" + ] + outer_channels = [ + AngularKetFJ(l_c=0, l_r=2, j_c=0.5, j_r=1.5, f_tot=1, species="Yb174"), + ] + + eigen_quantum_defects = [ + [0.75258093, 0.3826, -483.1], + ] + mixing_angles = [] + + +class Yb174_D2_HighN(FModel): + species_name = "Yb174_mqdt" + name = "D J=2, nu > 5" + f_tot = 2 + nu_range = (5.0, np.inf) + reference = "10.1103/PhysRevX.15.011009" + + inner_channels = [ + AngularKetLS(l_c=0, l_r=2, l_tot=2, s_tot=0, j_tot=2, species="Yb174"), # "6snd 1D2" + AngularKetLS(l_c=0, l_r=2, l_tot=2, s_tot=1, j_tot=2, species="Yb174"), # "6snd 3D2" + AngularKetDummy("4f13 5d 6snl a", f_tot=2), + AngularKetDummy("4f13 5d 6snl b", f_tot=2), + AngularKetLS(l_c=1, l_r=1, l_tot=2, s_tot=0, j_tot=2, species="Yb174"), # "6pnp 1D2" + ] + outer_channels = [ + AngularKetFJ(l_c=0, l_r=2, j_c=0.5, j_r=2.5, f_tot=2, species="Yb174"), + AngularKetFJ(l_c=0, l_r=2, j_c=0.5, j_r=1.5, f_tot=2, species="Yb174"), + AngularKetDummy("4f13 5d 6snl a", f_tot=2), + AngularKetDummy("4f13 5d 6snl b", f_tot=2), + AngularKetFJ(l_c=1, l_r=1, j_c=Unknown, j_r=Unknown, f_tot=2, species="Yb174"), # Jc could be either 1/2 or 3/2 + ] + manual_frame_transformation_outer_inner = np.array( + [ + [np.sqrt(3 / 5), np.sqrt(2 / 5), 0, 0, 0], + [-np.sqrt(2 / 5), np.sqrt(3 / 5), 0, 0, 0], + [0, 0, 1, 0, 0], + [0, 0, 0, 1, 0], + [0, 0, 0, 0, 1], + ] + ) + + eigen_quantum_defects = [ + [0.729513646, -0.0377841183], + [0.752292223, 0.104072325], + [0.19612036, 0], + [0.233752026, 0], + [0.152911249, 0], + ] + mixing_angles = [ + (0, 1, [0.21157531, -15.3844]), + (0, 2, 0.00521559431), + (0, 3, 0.0398131577), + (1, 3, -0.0071658109), + (0, 4, 0.10481227), + (1, 4, 0.0721660042), + ] + + +class Yb174_D3_HighN(FModel): + species_name = "Yb174_mqdt" + name = "D J=3, nu > 18" + f_tot = 3 + nu_range = (18.0, np.inf) + reference = "10.1103/PhysRevX.15.011009" + + inner_channels = [ + AngularKetLS(l_c=0, l_r=2, l_tot=2, s_tot=1, j_tot=3, species="Yb174"), # "6snd 3D3" + ] + outer_channels = [ + AngularKetFJ(l_c=0, l_r=2, j_c=0.5, j_r=2.5, f_tot=3, species="Yb174"), + ] + + eigen_quantum_defects = [ + [0.72902016, -0.705328923, 829.238844], + ] + mixing_angles = [] + + +class Yb174_F2_HighN(FModel): + species_name = "Yb174_mqdt" + name = "F J=2, nu > 25" + f_tot = 2 + nu_range = (25.0, np.inf) + reference = "arXiv:2507.11487v1" + + inner_channels = [ + AngularKetLS(l_c=0, l_r=3, l_tot=3, s_tot=1, j_tot=2, species="Yb174"), # "6snf 3F2" + ] + outer_channels = [ + AngularKetFJ(l_c=0, l_r=3, j_c=0.5, j_r=2.5, f_tot=2, species="Yb174"), + ] + + eigen_quantum_defects = [ + [0.0718252326, -1.00091963, -106.291066], + ] + mixing_angles = [] + + +class Yb174_F3_HighN(FModel): + species_name = "Yb174_mqdt" + name = "F J=3, nu > 7" + f_tot = 3 + nu_range = (7.0, np.inf) + reference = "10.1103/PhysRevA.112.042817" + + inner_channels = [ + AngularKetLS(l_c=0, l_r=3, l_tot=3, s_tot=0, j_tot=3, species="Yb174"), # "6snf 1F3" + AngularKetLS(l_c=0, l_r=3, l_tot=3, s_tot=1, j_tot=3, species="Yb174"), # "6snf 3F3" + AngularKetDummy("4f13 5d 6snl a", f_tot=3), + AngularKetDummy("4f13 5d 6snl b", f_tot=3), + AngularKetDummy("4f13 5d 6snl c", f_tot=3), + AngularKetDummy("4f13 5d 6snl d", f_tot=3), + AngularKetDummy("4f13 5d 6snl e", f_tot=3), + ] + outer_channels = [ + AngularKetFJ(l_c=0, l_r=3, j_c=0.5, j_r=3.5, f_tot=3, species="Yb174"), + AngularKetFJ(l_c=0, l_r=3, j_c=0.5, j_r=2.5, f_tot=3, species="Yb174"), + AngularKetDummy("4f13 5d 6snl a", f_tot=3), + AngularKetDummy("4f13 5d 6snl b", f_tot=3), + AngularKetDummy("4f13 5d 6snl c", f_tot=3), + AngularKetDummy("4f13 5d 6snl d", f_tot=3), + AngularKetDummy("4f13 5d 6snl e", f_tot=3), + ] + + eigen_quantum_defects = [ + [0.276158949, -12.7258012], + [0.0715123712, -0.768462937], + [0.239015576, 0], + [0.226770354, 0], + [0.175354845, 0], + [0.196660618, 0], + [0.21069642, 0], + ] + mixing_angles = [ + (0, 1, [0.0209955122, 0.251041249]), + (0, 2, -0.00411835457), + (0, 3, -0.0962784945), + (0, 4, 0.132826901), + (0, 5, -0.0439244317), + (0, 6, 0.0508460294), + (1, 2, -0.0376574252), + (1, 3, 0.026944623), + (1, 4, -0.0148474857), + (1, 5, -0.0521244126), + (1, 6, 0.0349516329), + ] + + +class Yb174_F4_HighN(FModel): + species_name = "Yb174_mqdt" + name = "F J=4, nu > 25" + f_tot = 4 + nu_range = (25.0, np.inf) + reference = "arXiv:2507.11487v1" + + inner_channels = [ + AngularKetLS(l_c=0, l_r=3, l_tot=3, s_tot=1, j_tot=4, species="Yb174"), # "6snf 3F4" + ] + outer_channels = [ + AngularKetFJ(l_c=0, l_r=3, j_c=0.5, j_r=3.5, f_tot=4, species="Yb174"), + ] + + eigen_quantum_defects = [ + [0.0839027969, -2.91009023], + ] + mixing_angles = [] + + +class Yb174_G3_HighN(FModel): + species_name = "Yb174_mqdt" + name = "G J=3, nu > 25" + f_tot = 3 + nu_range = (25.0, np.inf) + reference = "arXiv:2507.11487v1" + + inner_channels = [ + AngularKetLS(l_c=0, l_r=4, l_tot=4, s_tot=1, j_tot=3, species="Yb174"), # "6sng 3G3" + ] + outer_channels = [ + AngularKetFJ(l_c=0, l_r=4, j_c=0.5, j_r=3.5, f_tot=3, species="Yb174"), + ] + + eigen_quantum_defects = [ + [0.0260964574, -0.14139526], + ] + mixing_angles = [] + + +class Yb174_G4_HighN(FModel): + species_name = "Yb174_mqdt" + name = "G J=4, nu > 25" + f_tot = 4 + nu_range = (25.0, np.inf) + reference = "10.1103/PhysRevA.112.042817" + + inner_channels = [ + AngularKetJJ(l_c=0, l_r=4, j_c=0.5, j_r=4.5, j_tot=4, species="Yb174"), # "6sng +G4" + AngularKetJJ(l_c=0, l_r=4, j_c=0.5, j_r=3.5, j_tot=4, species="Yb174"), # "6sng -G4" + ] + outer_channels = [ + AngularKetFJ(l_c=0, l_r=4, j_c=0.5, j_r=4.5, f_tot=4, species="Yb174"), + AngularKetFJ(l_c=0, l_r=4, j_c=0.5, j_r=3.5, f_tot=4, species="Yb174"), + ] + + eigen_quantum_defects = [ + [0.0262659964, -0.148808463], + [0.0254568575, -0.134219071], + ] + mixing_angles = [ + (0, 1, -0.089123698), + ] + + +class Yb174_G5_HighN(FModel): + species_name = "Yb174_mqdt" + name = "G J=5, nu > 25" + f_tot = 5 + nu_range = (25.0, np.inf) + reference = "10.1103/PhysRevA.112.042817" + + inner_channels = [ + AngularKetLS(l_c=0, l_r=4, l_tot=4, s_tot=1, j_tot=5, species="Yb174"), # "6sng 3G5" + ] + outer_channels = [ + AngularKetFJ(l_c=0, l_r=4, j_c=0.5, j_r=4.5, f_tot=5, species="Yb174"), + ] + + eigen_quantum_defects = [ + [0.02536571, -0.18507079], + ] + mixing_angles = [] + + +# -------------------------------------------------------- +# MQDT models valid at small n +# -------------------------------------------------------- + + +class Yb174_S0_LowN(FModel): + species_name = "Yb174_mqdt" + name = "S J=0, 1 < nu < 2" + f_tot = 0 + nu_range = (1.0, 2.0) + reference = None + + inner_channels = [ + AngularKetLS(l_c=0, l_r=0, l_tot=0, s_tot=0, j_tot=0, species="Yb174"), # "6sns 1S0" + ] + outer_channels = [ + AngularKetFJ(l_c=0, l_r=0, j_c=0.5, j_r=0.5, f_tot=0, species="Yb174"), + ] + + eigen_quantum_defects = [ + [0.525055, 0], + ] + mixing_angles = [] + + +class Yb174_S1_LowN(FModel): + species_name = "Yb174_mqdt" + name = "S J=1, 2 < nu < 26" + f_tot = 1 + nu_range = (2.0, 26.0) + reference = None + + inner_channels = [ + AngularKetLS(l_c=0, l_r=0, l_tot=0, s_tot=1, j_tot=1, species="Yb174"), # "6sns 3S1" + ] + outer_channels = [ + AngularKetFJ(l_c=0, l_r=0, j_c=0.5, j_r=0.5, f_tot=1, species="Yb174"), + ] + + eigen_quantum_defects = [ + [0.432841, 0.724559, -1.95424], + ] + mixing_angles = [] + + +class Yb174_P0_LowN(FModel): + species_name = "Yb174_mqdt" + name = "P J=0, 1.5 < nu < 5.5" + f_tot = 0 + nu_range = (1.5, 5.5) + reference = None + + inner_channels = [ + AngularKetLS(l_c=0, l_r=1, l_tot=1, s_tot=1, j_tot=0, species="Yb174"), # "6snp 3P0" + ] + outer_channels = [ + AngularKetFJ(l_c=0, l_r=1, j_c=0.5, j_r=0.5, f_tot=0, species="Yb174"), + ] + + eigen_quantum_defects = [ + [0.969279, 0.288219, 1.36228], + ] + mixing_angles = [] + + +class Yb174_P1_Lowest(FModel): + species_name = "Yb174_mqdt" + name = "P J=1, 1.7 < nu < 2.7" + f_tot = 1 + nu_range = (1.7, 2.7) + reference = None + + inner_channels = [ + AngularKetLS(l_c=0, l_r=1, l_tot=1, s_tot=0, j_tot=1, species="Yb174"), # "6snp 1P1" + AngularKetLS(l_c=0, l_r=1, l_tot=1, s_tot=1, j_tot=1, species="Yb174"), # "6snp 3P1" + ] + outer_channels = [ + AngularKetFJ(l_c=0, l_r=1, j_c=0.5, j_r=1.5, f_tot=1, species="Yb174"), + AngularKetFJ(l_c=0, l_r=1, j_c=0.5, j_r=0.5, f_tot=1, species="Yb174"), + ] + + eigen_quantum_defects = [ + [0.161083, 0], + [0.920424, 0], + ] + mixing_angles = [ + (0, 1, [-0.426128, 6.272986]), + ] + + +class Yb174_P1_LowN(FModel): + species_name = "Yb174_mqdt" + name = "P J=1, 2.7 < nu < 5.7" + f_tot = 1 + nu_range = (2.7, 5.7) + reference = None + + inner_channels = [ + AngularKetLS(l_c=0, l_r=1, l_tot=1, s_tot=0, j_tot=1, species="Yb174"), # "6snp 1P1" + AngularKetLS(l_c=0, l_r=1, l_tot=1, s_tot=1, j_tot=1, species="Yb174"), # "6snp 3P1" + ] + outer_channels = [ + AngularKetFJ(l_c=0, l_r=1, j_c=0.5, j_r=1.5, f_tot=1, species="Yb174"), + AngularKetFJ(l_c=0, l_r=1, j_c=0.5, j_r=0.5, f_tot=1, species="Yb174"), + ] + + eigen_quantum_defects = [ + [0.967223, -3.03997, 0.569205], + [0.967918, 0.25116, 0.868505], + ] + mixing_angles = [] + + +class Yb174_P2_LowN(FModel): + species_name = "Yb174_mqdt" + name = "P J=2, 1.5 < nu < 4.5" + f_tot = 2 + nu_range = (1.5, 4.5) + reference = None + + inner_channels = [ + AngularKetLS(l_c=0, l_r=1, l_tot=1, s_tot=1, j_tot=2, species="Yb174"), # "6snp 3P2" + ] + outer_channels = [ + AngularKetFJ(l_c=0, l_r=1, j_c=0.5, j_r=1.5, f_tot=2, species="Yb174"), + ] + + eigen_quantum_defects = [ + [0.906105, 0.383471, 1.23512], + ] + mixing_angles = [] + + +class Yb174_D1_LowN(FModel): + species_name = "Yb174_mqdt" + name = "D J=1, 2 < nu < 26" + f_tot = 1 + nu_range = (2.0, 26.0) + reference = None + + inner_channels = [ + AngularKetLS(l_c=0, l_r=2, l_tot=2, s_tot=1, j_tot=1, species="Yb174"), # "6snd 3D1" + ] + outer_channels = [ + AngularKetFJ(l_c=0, l_r=2, j_c=0.5, j_r=1.5, f_tot=1, species="Yb174"), + ] + + eigen_quantum_defects = [ + [0.758222, -0.017906, 3.392161], + ] + mixing_angles = [] + + +class Yb174_D2_LowN(FModel): + species_name = "Yb174_mqdt" + name = "D J=2, 2 < nu < 5" + f_tot = 2 + nu_range = (2.0, 5.0) + reference = None + + inner_channels = [ + AngularKetLS(l_c=0, l_r=2, l_tot=2, s_tot=0, j_tot=2, species="Yb174"), # "6snd 1D2" + AngularKetLS(l_c=0, l_r=2, l_tot=2, s_tot=1, j_tot=2, species="Yb174"), # "6snd 3D2" + ] + outer_channels = [ + AngularKetFJ(l_c=0, l_r=2, j_c=0.5, j_r=2.5, f_tot=2, species="Yb174"), + AngularKetFJ(l_c=0, l_r=2, j_c=0.5, j_r=1.5, f_tot=2, species="Yb174"), + ] + + eigen_quantum_defects = [ + [0.703156, 0.973192], + [0.724546, 0.372621], + ] + mixing_angles = [ + (0, 1, [0.948409, 2.121270]), + ] + + +class Yb174_D3_LowN(FModel): + species_name = "Yb174_mqdt" + name = "D J=3, 2 < nu < 18" + f_tot = 3 + nu_range = (2.0, 18.0) + reference = None + + inner_channels = [ + AngularKetLS(l_c=0, l_r=2, l_tot=2, s_tot=1, j_tot=3, species="Yb174"), # "6snd 3D3" + ] + outer_channels = [ + AngularKetFJ(l_c=0, l_r=2, j_c=0.5, j_r=2.5, f_tot=3, species="Yb174"), + ] + + eigen_quantum_defects = [ + [0.734512, -0.019501, 3.459114], + ] + mixing_angles = [] From 2fecdc5e72e7abf5fd1fb9abbef60fdd895fc72d Mon Sep 17 00:00:00 2001 From: johannes-moegerle Date: Tue, 10 Mar 2026 12:09:30 +0100 Subject: [PATCH 05/18] TODO fixup species mqdt --- src/rydstate/angular/angular_ket.py | 4 ++-- 1 file changed, 2 insertions(+), 2 deletions(-) diff --git a/src/rydstate/angular/angular_ket.py b/src/rydstate/angular/angular_ket.py index 74ecb58..3b302a3 100644 --- a/src/rydstate/angular/angular_ket.py +++ b/src/rydstate/angular/angular_ket.py @@ -95,9 +95,9 @@ def __init__( """ if species is not None: if isinstance(species, str): - from rydstate.species import SpeciesObjectSQDT # noqa: PLC0415 + from rydstate.species.sqdt import SpeciesObjectSQDT # noqa: PLC0415 - species = SpeciesObjectSQDT.from_name(species) + species = SpeciesObjectSQDT.from_name(species.replace("_mqdt", "")) # use i_c = 0 for species without defined nuclear spin (-> ignore hyperfine) if i_c is not None and i_c != species.i_c_number: raise ValueError(f"Nuclear spin i_c={i_c} does not match the species {species} with i_c={species.i_c}.") From 9c7deb9b2f508123db6628f2927b69dec7aa3cdb Mon Sep 17 00:00:00 2001 From: johannes-moegerle Date: Tue, 10 Mar 2026 12:13:17 +0100 Subject: [PATCH 06/18] [tests] add tests for fmodel --- tests/test_mqdt_models.py | 159 ++++++++++++++++++++++++++++++++++++++ 1 file changed, 159 insertions(+) create mode 100644 tests/test_mqdt_models.py diff --git a/tests/test_mqdt_models.py b/tests/test_mqdt_models.py new file mode 100644 index 0000000..56b98d1 --- /dev/null +++ b/tests/test_mqdt_models.py @@ -0,0 +1,159 @@ +from __future__ import annotations + +import re + +import numpy as np +import pytest +import rydstate.species.mqdt.sr87 +import rydstate.species.mqdt.sr88 +import rydstate.species.mqdt.yb171 +import rydstate.species.mqdt.yb173 +import rydstate.species.mqdt.yb174 # noqa: F401 +from rydstate.angular.utils import is_dummy_ket +from rydstate.species.mqdt.fmodel import FModel + + +def _all_fmodels() -> list[FModel]: + """Collect all concrete FModel subclasses.""" + return [cls() for cls in FModel.__subclasses__() if hasattr(cls, "name") and cls.name is not None] + + +ALL_MODELS = _all_fmodels() + + +@pytest.fixture(params=ALL_MODELS, ids=lambda cls: cls.name) +def model(request: pytest.FixtureRequest) -> FModel: + return request.param # type: ignore[no-any-return] + + +def test_all_models_discovered() -> None: + """Sanity check: we should find at least 80 FModel subclasses.""" + assert len(ALL_MODELS) >= 80 + + +def test_channel_count_consistency(model: FModel) -> None: + """Inner channels, outer channels, and eigen quantum defects must have the same length.""" + n_inner = len(model.inner_channels) + n_outer = len(model.outer_channels) + n_defects = len(model.eigen_quantum_defects) + assert n_inner == n_outer, f"{model.name}: inner_channels ({n_inner}) != outer_channels ({n_outer})" + assert n_inner == n_defects, f"{model.name}: channels ({n_inner}) != eigen_quantum_defects ({n_defects})" + + +def test_nu_range_valid(model: FModel) -> None: + """nu_range must be a 2-tuple with min < max.""" + nu_min, nu_max = model.nu_range + assert nu_min < nu_max, f"{model.name}: nu_range {model.nu_range} has min >= max" + assert nu_min > 0, f"{model.name}: nu_range min must be positive" + + +def test_f_tot_consistency(model: FModel) -> None: + """All non-dummy channels must have f_tot matching the model's f_tot.""" + for i, ch in enumerate(model.inner_channels): + if not is_dummy_ket(ch): + assert ch.f_tot == model.f_tot, ( + f"{model.name}: inner_channels[{i}].f_tot={ch.f_tot} != model.f_tot={model.f_tot}" + ) + for i, och in enumerate(model.outer_channels): + if not is_dummy_ket(och): + assert och.f_tot == model.f_tot, ( + f"{model.name}: outer_channels[{i}].f_tot={och.f_tot} != model.f_tot={model.f_tot}" + ) + + +def test_parity_consistency(model: FModel) -> None: + """All non-dummy channels must have the same parity (-1)^(l_c+l_r).""" + parities: list[int] = [(-1) ** (ch.l_c + ch.l_r) for ch in model.inner_channels if not is_dummy_ket(ch)] + parities.extend((-1) ** (och.l_c + och.l_r) for och in model.outer_channels if not is_dummy_ket(och)) + assert len(set(parities)) <= 1, f"{model.name}: channels have inconsistent parity" + + +def test_dummy_channels_match(model: FModel) -> None: + """Dummy channels must appear at the same positions in inner and outer channel lists.""" + for i, (inner, outer) in enumerate(zip(model.inner_channels, model.outer_channels, strict=True)): + inner_is_dummy = is_dummy_ket(inner) + outer_is_dummy = is_dummy_ket(outer) + assert inner_is_dummy == outer_is_dummy, ( + f"{model.name}: channel {i} dummy mismatch " + f"(inner={'dummy' if inner_is_dummy else 'real'}, outer={'dummy' if outer_is_dummy else 'real'})" + ) + + +def test_mixing_angles_indices_valid(model: FModel) -> None: + """Mixing angle indices must refer to valid channel positions.""" + n_channels = len(model.inner_channels) + for entry in model.mixing_angles: + i, j = entry[0], entry[1] + assert 0 <= i < n_channels, f"{model.name}: mixing_angles index {i} out of range [0, {n_channels})" + assert 0 <= j < n_channels, f"{model.name}: mixing_angles index {j} out of range [0, {n_channels})" + assert i != j, f"{model.name}: mixing_angles has self-coupling ({i}, {j})" + + +def test_model_name_contains_quantum_number(model: FModel) -> None: + """Model name must contain F=X/Y or J=X matching the model's f_tot.""" + assert model.name is not None + # Match F=X/Y or J=X patterns (integers and fractions) + match = re.search(r"[FJ]=(\d+/?\.?\d*)", model.name) + assert match is not None, f"{model.name}: name '{model.name}' does not contain F=... or J=..." + s = match.group(1) + if "/" in s: + num, den = s.split("/") + f_val = float(num) / float(den) + else: + f_val = float(s) + assert f_val == model.f_tot, f"{model.name}: name says F/J={f_val} but f_tot={model.f_tot}" + + +def test_reference_field_set(model: FModel) -> None: + """Every model must have an explicit reference field (str or None).""" + assert hasattr(model, "reference"), f"{model.name}: missing 'reference' field" + + +def test_species_field_set(model: FModel) -> None: + """Every model must have a species field.""" + assert model.species_name is not None, f"{model.name}: species is None" + + +def test_eigen_quantum_defects_format(model: FModel) -> None: + """Each eigen quantum defect entry must be a list/tuple of numeric values.""" + for i, defect in enumerate(model.eigen_quantum_defects): + if isinstance(defect, (list, tuple)): + for j, val in enumerate(defect): + assert isinstance(val, (int, float)), ( + f"{model.name}: eigen_quantum_defects[{i}][{j}] is {type(val)}, expected numeric" + ) + else: + assert isinstance(defect, (int, float)), ( + f"{model.name}: eigen_quantum_defects[{i}] is {type(defect)}, expected numeric or list" + ) + + +def test_inner_outer_channel_l_r_match(model: FModel) -> None: + """Non-dummy inner and outer channels at the same position must have the same l_r.""" + for i, (inner, outer) in enumerate(zip(model.inner_channels, model.outer_channels, strict=True)): + if is_dummy_ket(inner) or is_dummy_ket(outer): + continue + assert inner.l_r == outer.l_r, ( + f"{model.name}: channel {i} l_r mismatch (inner.l_r={inner.l_r}, outer.l_r={outer.l_r})" + ) + + +def test_at_least_one_real_channel(model: FModel) -> None: + """Every model must have at least one non-dummy channel.""" + real_channels = [ch for ch in model.inner_channels if not is_dummy_ket(ch)] + assert len(real_channels) >= 1, f"{model.name}: no real (non-dummy) channels" + + +def test_inner_outer_unitary(model: FModel) -> None: + """The frame transformation matrix from inner to outer channels must be unitary.""" + unitary = model.calc_frame_transformation_outer_inner() + msg = f"{model.species_name} - {model.name}: frame transformation (outer - inner) is not unitary" + np.testing.assert_allclose(unitary.conj().T @ unitary, np.eye(unitary.shape[0]), atol=1e-10, err_msg=msg) + + rotation = model.calc_frame_transformation_inner_closecoupling(nu_ref=30.5) + msg = f"{model.species_name} - {model.name}: frame transformation (inner - closecoupling) is not unitary" + np.testing.assert_allclose(rotation.conj().T @ rotation, np.eye(rotation.shape[0]), atol=1e-10, err_msg=msg) + + full = model.calc_frame_transformation(nu_ref=30.5) + msg = f"{model.species_name} - {model.name}: full frame transformation U=QR is not unitary" + np.testing.assert_allclose(full.conj().T @ full, np.eye(full.shape[0]), atol=1e-10, err_msg=msg) From eb9620c668fd949ea15ba2a699b1dfa47cf69c13 Mon Sep 17 00:00:00 2001 From: johannes-moegerle Date: Tue, 10 Mar 2026 12:11:30 +0100 Subject: [PATCH 07/18] allow mqdt species in RydbergStateSQDT --- src/rydstate/rydberg/rydberg_sqdt.py | 18 +++++++++++++----- 1 file changed, 13 insertions(+), 5 deletions(-) diff --git a/src/rydstate/rydberg/rydberg_sqdt.py b/src/rydstate/rydberg/rydberg_sqdt.py index a6934ac..050a8af 100644 --- a/src/rydstate/rydberg/rydberg_sqdt.py +++ b/src/rydstate/rydberg/rydberg_sqdt.py @@ -13,6 +13,7 @@ from rydstate.radial import RadialKet from rydstate.rydberg.rydberg_base import RydbergStateBase from rydstate.species import SpeciesObjectSQDT +from rydstate.species.species_object import SpeciesObject from rydstate.species.utils import calc_energy_from_nu from rydstate.units import BaseQuantities, MatrixElementOperatorRanks, ureg @@ -27,7 +28,7 @@ class RydbergStateSQDT(RydbergStateBase): - species: SpeciesObjectSQDT + species: SpeciesObject """The atomic species of the Rydberg state.""" angular: AngularKetBase @@ -35,7 +36,7 @@ class RydbergStateSQDT(RydbergStateBase): def __init__( self, - species: str | SpeciesObjectSQDT, + species: str | SpeciesObject, n: int | None = None, nu: float | None = None, s_c: float | None = None, @@ -74,7 +75,7 @@ def __init__( """ if isinstance(species, str): - species = SpeciesObjectSQDT.from_name(species) + species = SpeciesObject.from_name(species) self.species = species self.angular = quantum_numbers_to_angular_ket( @@ -106,7 +107,7 @@ def _set_qn_as_attributes(self) -> None: @classmethod def from_angular_ket( cls: type[Self], - species: str | SpeciesObjectSQDT, + species: str | SpeciesObject, angular_ket: AngularKetBase, n: int | None = None, nu: float | None = None, @@ -115,7 +116,7 @@ def from_angular_ket( obj = cls.__new__(cls) if isinstance(species, str): - species = SpeciesObjectSQDT.from_name(species) + species = SpeciesObject.from_name(species) obj.species = species obj.n = n @@ -425,6 +426,9 @@ def _get_transition_rates_au( raise NotImplementedError( "For alkaline earth atoms transition rates are only implemented for LS coupling scheme." ) + if not isinstance(self.species, SpeciesObjectSQDT): + raise NotImplementedError("transition rates are currently only implemented for SQDT species.") + from rydstate.basis import BasisSQDTAlkali, BasisSQDTAlkalineLS # noqa: PLC0415 basis_class = BasisSQDTAlkali if self.species.number_valence_electrons == 1 else BasisSQDTAlkalineLS @@ -533,6 +537,7 @@ class RydbergStateSQDTAlkali(RydbergStateSQDT): """Create an Alkali Rydberg state, including the radial and angular states.""" angular: AngularKetLS + species: SpeciesObjectSQDT def __init__( self, @@ -580,6 +585,7 @@ class RydbergStateSQDTAlkalineLS(RydbergStateSQDT): """Create an Alkaline Rydberg state, including the radial and angular states.""" angular: AngularKetLS + species: SpeciesObjectSQDT def __init__( self, @@ -629,6 +635,7 @@ class RydbergStateSQDTAlkalineJJ(RydbergStateSQDT): """Create an Alkaline Rydberg state, including the radial and angular states.""" angular: AngularKetJJ + species: SpeciesObjectSQDT def __init__( self, @@ -678,6 +685,7 @@ class RydbergStateSQDTAlkalineFJ(RydbergStateSQDT): """Create an Alkaline Rydberg state, including the radial and angular states.""" angular: AngularKetFJ + species: SpeciesObjectSQDT def __init__( self, From 08e735154a1c1e57d86f13b52101ad3c2f0607ba Mon Sep 17 00:00:00 2001 From: johannes-moegerle Date: Tue, 10 Mar 2026 12:11:42 +0100 Subject: [PATCH 08/18] add rydberg state dummy --- src/rydstate/rydberg/rydberg_dummy.py | 26 ++++++++++++++++++++++++++ 1 file changed, 26 insertions(+) create mode 100644 src/rydstate/rydberg/rydberg_dummy.py diff --git a/src/rydstate/rydberg/rydberg_dummy.py b/src/rydstate/rydberg/rydberg_dummy.py new file mode 100644 index 0000000..6d89360 --- /dev/null +++ b/src/rydstate/rydberg/rydberg_dummy.py @@ -0,0 +1,26 @@ +from __future__ import annotations + +from typing import TYPE_CHECKING, TypeGuard + +from rydstate.rydberg.rydberg_sqdt import RydbergStateSQDT + +if TYPE_CHECKING: + from rydstate.angular.angular_ket_dummy import AngularKetDummy + + +class RydbergStateSQDTDummy(RydbergStateSQDT): + """Create an Alkaline Rydberg state, including the radial and angular states.""" + + angular: AngularKetDummy + + def __init__(self) -> None: + raise RuntimeError("RydbergStateSQDTDummy ican only be created via the from_angular_ket method.") + + def __repr__(self) -> str: + species, nu = self.species, self.nu + return f"{self.__class__.__name__}({species.name}, {nu=}, angular_ket={self.angular})" + + +def is_dummy_rydberg_state(state: RydbergStateSQDT) -> TypeGuard[RydbergStateSQDTDummy]: + """Check if state is a RydbergStateSQDTDummy.""" + return isinstance(state, RydbergStateSQDTDummy) From 4c7489e600d625f824cc7aaebd3c3c590414e2ce Mon Sep 17 00:00:00 2001 From: johannes-moegerle Date: Tue, 10 Mar 2026 12:12:30 +0100 Subject: [PATCH 09/18] add mqdt stuff (state, basis, calculations) --- src/rydstate/__init__.py | 10 +- src/rydstate/basis/__init__.py | 10 +- src/rydstate/basis/basis_mqdt.py | 73 +++++++++++ src/rydstate/rydberg/__init__.py | 2 + src/rydstate/rydberg/rydberg_mqdt.py | 184 +++++++++++++++++++++++++++ src/rydstate/utils/__init__.py | 0 src/rydstate/utils/linalg.py | 116 +++++++++++++++++ 7 files changed, 388 insertions(+), 7 deletions(-) create mode 100644 src/rydstate/basis/basis_mqdt.py create mode 100644 src/rydstate/rydberg/rydberg_mqdt.py create mode 100644 src/rydstate/utils/__init__.py create mode 100644 src/rydstate/utils/linalg.py diff --git a/src/rydstate/__init__.py b/src/rydstate/__init__.py index 739e3ad..42ffe79 100644 --- a/src/rydstate/__init__.py +++ b/src/rydstate/__init__.py @@ -1,11 +1,7 @@ from rydstate import angular, basis, radial, rydberg, species -from rydstate.basis import ( - BasisSQDTAlkali, - BasisSQDTAlkalineFJ, - BasisSQDTAlkalineJJ, - BasisSQDTAlkalineLS, -) +from rydstate.basis import BasisMQDT, BasisSQDTAlkali, BasisSQDTAlkalineFJ, BasisSQDTAlkalineJJ, BasisSQDTAlkalineLS from rydstate.rydberg import ( + RydbergStateMQDT, RydbergStateSQDT, RydbergStateSQDTAlkali, RydbergStateSQDTAlkalineFJ, @@ -15,10 +11,12 @@ from rydstate.units import ureg __all__ = [ + "BasisMQDT", "BasisSQDTAlkali", "BasisSQDTAlkalineFJ", "BasisSQDTAlkalineJJ", "BasisSQDTAlkalineLS", + "RydbergStateMQDT", "RydbergStateSQDT", "RydbergStateSQDTAlkali", "RydbergStateSQDTAlkalineFJ", diff --git a/src/rydstate/basis/__init__.py b/src/rydstate/basis/__init__.py index 7f15364..93ef96c 100644 --- a/src/rydstate/basis/__init__.py +++ b/src/rydstate/basis/__init__.py @@ -1,4 +1,12 @@ from rydstate.basis.basis_base import BasisBase +from rydstate.basis.basis_mqdt import BasisMQDT from rydstate.basis.basis_sqdt import BasisSQDTAlkali, BasisSQDTAlkalineFJ, BasisSQDTAlkalineJJ, BasisSQDTAlkalineLS -__all__ = ["BasisBase", "BasisSQDTAlkali", "BasisSQDTAlkalineFJ", "BasisSQDTAlkalineJJ", "BasisSQDTAlkalineLS"] +__all__ = [ + "BasisBase", + "BasisMQDT", + "BasisSQDTAlkali", + "BasisSQDTAlkalineFJ", + "BasisSQDTAlkalineJJ", + "BasisSQDTAlkalineLS", +] diff --git a/src/rydstate/basis/basis_mqdt.py b/src/rydstate/basis/basis_mqdt.py new file mode 100644 index 0000000..c6726c2 --- /dev/null +++ b/src/rydstate/basis/basis_mqdt.py @@ -0,0 +1,73 @@ +from __future__ import annotations + +import logging + +import numpy as np + +from rydstate.angular import AngularKetFJ +from rydstate.basis.basis_base import BasisBase +from rydstate.rydberg import RydbergStateMQDT +from rydstate.rydberg.rydberg_mqdt import get_mqdt_states_from_fmodel +from rydstate.species import FModelSQDT, SpeciesObjectMQDT + +logger = logging.getLogger(__name__) + + +class BasisMQDT(BasisBase[RydbergStateMQDT]): + def __init__( # noqa: C901, PLR0912 + self, + species: str | SpeciesObjectMQDT, + nu_min: int = 0, + nu_max: int | None = None, + *, + f_tot: float | tuple[float, float] | None = None, + skip_high_l: bool = True, + n_min_high_l: int = 0, + ) -> None: + if isinstance(species, str): + species = SpeciesObjectMQDT.from_name(species) + self.species = species + + if nu_max is None: + raise ValueError("nu_max must be given") + + self.models = [] + i_c, j_c, s_r = self.species.i_c, 0.5, 0.5 + self.states = [] + for l_r in range(int(nu_max) + 1): + for j_r in np.arange(abs(l_r - s_r), l_r + s_r + 1): + for f_c in np.arange(abs(j_c - i_c), j_c + i_c + 1): + for _f_tot in np.arange(abs(f_c - j_r), f_c + j_r + 1): + if f_tot is not None: + if isinstance(f_tot, tuple): + if not (f_tot[0] <= _f_tot <= f_tot[1]): + continue + elif _f_tot != f_tot: + continue + channel = AngularKetFJ( + l_r=l_r, j_r=float(j_r), f_c=float(f_c), f_tot=float(_f_tot), species=species + ) + self.models.extend(species.get_mqdt_models(channel)) + self.models = list(set(self.models)) # remove duplicates + + logger.debug("Calculating MQDT states...") + self.states: list[RydbergStateMQDT] = [] + for model in self.models: + _nu_min = nu_min + if isinstance(model, FModelSQDT): + if skip_high_l: + continue + _nu_min = n_min_high_l + logger.debug(" calculating states for model %s with nu_min=%s, nu_max=%s", model.name, _nu_min, nu_max) + _states = get_mqdt_states_from_fmodel(model, _nu_min, nu_max) + if len(_states) == 0: + logger.debug(" no states found for model %s", model.name) + else: + logger.debug( + " model %s: nu_min=%s, nu_max=%s, total states=%d", + model.name, + _states[0].nu_ref, + _states[-1].nu_ref, + len(_states), + ) + self.states.extend(_states) diff --git a/src/rydstate/rydberg/__init__.py b/src/rydstate/rydberg/__init__.py index f7a895b..b13d7a7 100644 --- a/src/rydstate/rydberg/__init__.py +++ b/src/rydstate/rydberg/__init__.py @@ -1,4 +1,5 @@ from rydstate.rydberg.rydberg_base import RydbergStateBase +from rydstate.rydberg.rydberg_mqdt import RydbergStateMQDT from rydstate.rydberg.rydberg_sqdt import ( RydbergStateSQDT, RydbergStateSQDTAlkali, @@ -9,6 +10,7 @@ __all__ = [ "RydbergStateBase", + "RydbergStateMQDT", "RydbergStateSQDT", "RydbergStateSQDTAlkali", "RydbergStateSQDTAlkalineFJ", diff --git a/src/rydstate/rydberg/rydberg_mqdt.py b/src/rydstate/rydberg/rydberg_mqdt.py new file mode 100644 index 0000000..6680267 --- /dev/null +++ b/src/rydstate/rydberg/rydberg_mqdt.py @@ -0,0 +1,184 @@ +from __future__ import annotations + +import logging +from typing import TYPE_CHECKING, overload + +import numpy as np + +from rydstate.angular import AngularState +from rydstate.angular.utils import is_dummy_ket +from rydstate.rydberg.rydberg_base import RydbergStateBase +from rydstate.rydberg.rydberg_dummy import RydbergStateSQDTDummy, is_dummy_rydberg_state +from rydstate.rydberg.rydberg_sqdt import RydbergStateSQDT, RydbergStateSQDTAlkalineFJ +from rydstate.species import FModel +from rydstate.utils.linalg import calc_nullvector + +if TYPE_CHECKING: + from collections.abc import Iterator, Sequence + + from rydstate.angular import AngularKetFJ + from rydstate.angular.angular_ket_dummy import AngularKetDummy + from rydstate.species import FModel + from rydstate.units import MatrixElementOperator, NDArray, PintFloat + + +logger = logging.getLogger(__name__) + + +class RydbergStateMQDT(RydbergStateBase): + angular: AngularState[AngularKetFJ | AngularKetDummy] + """Return the angular part of the MQDT state as an AngularState.""" + + def __init__( + self, + coefficients: Sequence[float] | NDArray, + sqdt_states: Sequence[RydbergStateSQDTAlkalineFJ | RydbergStateSQDTDummy], + nu_ref: float, + *, + warn_if_not_normalized: bool = True, + ) -> None: + self.coefficients = np.array(coefficients) + self.sqdt_states = sqdt_states + self.species = sqdt_states[0].species + self._nu_ref = nu_ref + self.angular = AngularState(self.coefficients.tolist(), [ket.angular for ket in sqdt_states]) + + if len(coefficients) != len(sqdt_states): + raise ValueError("Length of coefficients and sqdt_states must be the same.") + if not all( + (type(sqdt_state) is type(sqdt_states[0])) or is_dummy_rydberg_state(sqdt_state) + for sqdt_state in sqdt_states + ): + raise ValueError("All sqdt_states must be of the same type.") + if not all((sqdt_state.species is sqdt_states[0].species) for sqdt_state in sqdt_states): + raise ValueError("All sqdt_states must be of the same species.") + if len(set(sqdt_states)) != len(sqdt_states): + raise ValueError("RydbergStateMQDT initialized with duplicate sqdt_states.") + if abs(self.norm - 1) > 1e-10 and warn_if_not_normalized: + logger.warning( + "RydbergStateMQDT initialized with non-normalized coefficients " + "(norm=%s, coefficients=%s, sqdt_states=%s)", + self.norm, + coefficients, + sqdt_states, + ) + if self.norm > 1: + self.coefficients /= self.norm + + def __iter__(self) -> Iterator[tuple[float, RydbergStateSQDTAlkalineFJ | RydbergStateSQDTDummy]]: + return zip(self.coefficients, self.sqdt_states, strict=True).__iter__() + + def __repr__(self) -> str: + terms = [f"{coeff}*{sqdt_state!r}" for coeff, sqdt_state in self] + return f"{self.__class__.__name__}({', '.join(terms)})" + + def __str__(self) -> str: + terms = [f"{coeff}*{sqdt_state!s}" for coeff, sqdt_state in self] + return f"{', '.join(terms)}" + + @property + def nu_ref(self) -> float: + return self._nu_ref + + @property + def norm(self) -> float: + """Return the norm of the state (should be 1).""" + return float(np.linalg.norm(self.coefficients)) + + def calc_reduced_overlap(self, other: RydbergStateBase) -> float: + """Calculate the reduced overlap (ignoring the magnetic quantum number m).""" + other_iter: list[tuple[float, RydbergStateSQDT]] + if isinstance(other, RydbergStateSQDT): + other_iter = [(1.0, other)] + elif isinstance(other, RydbergStateMQDT): + other_iter = [(coeff, sqdt) for coeff, sqdt in other] + else: + raise NotImplementedError(f"calc_reduced_overlap not implemented for {type(self)=}, {type(other)=}") + + ov = 0 + for coeff1, sqdt1 in self: + for coeff2, sqdt2 in other_iter: + ov += np.conjugate(coeff1) * coeff2 * sqdt1.calc_reduced_overlap(sqdt2) + return ov + + @overload # type: ignore [override] + def calc_reduced_matrix_element( + self, other: RydbergStateBase, operator: MatrixElementOperator, unit: None = None + ) -> PintFloat: ... + + @overload + def calc_reduced_matrix_element( + self, other: RydbergStateBase, operator: MatrixElementOperator, unit: str + ) -> float: ... + + def calc_reduced_matrix_element( + self, other: RydbergStateBase, operator: MatrixElementOperator, unit: str | None = None + ) -> PintFloat | float: + r"""Calculate the reduced angular matrix element. + + This means, calculate the following matrix element: + + .. math:: + \left\langle self || \hat{O}^{(\kappa)} || other \right\rangle + + """ + other_iter: list[tuple[float, RydbergStateSQDT]] + if isinstance(other, RydbergStateSQDT): + other_iter = [(1.0, other)] + elif isinstance(other, RydbergStateMQDT): + other_iter = [(coeff, sqdt) for coeff, sqdt in other] + else: + raise NotImplementedError(f"calc_reduced_matrix_element not implemented for {type(self)=}, {type(other)=}") + + value = 0 + for coeff1, sqdt1 in self: + for coeff2, sqdt2 in other_iter: + value += np.conjugate(coeff1) * coeff2 * sqdt1.calc_reduced_matrix_element(sqdt2, operator, unit=unit) + return value + + +def get_mqdt_states_from_fmodel( + model: FModel, + nu_min: float | None = None, + nu_max: float | None = None, + *, + overwrite_model_limits: bool = False, +) -> list[RydbergStateMQDT]: + if overwrite_model_limits: + if nu_min is None or nu_max is None: + raise ValueError("nu_min and nu_max must be given if overwrite_model_limits is True.") + else: + nu_min = model.nu_min if nu_min is None else max(nu_min, model.nu_min) + nu_max = model.nu_max if nu_max is None else min(nu_max, model.nu_max) + if np.isinf(nu_max): + raise ValueError("nu_max must be finite to calculate MQDT states.") + + nu_ref_list = model.calc_detm_roots(nu_min, nu_max) + if len(nu_ref_list) == 0: + logger.warning( + "No MQDT states found in the range nu_min=%s, nu_max=%s for model %s", nu_min, nu_max, model.name + ) + return [] + + states: list[RydbergStateMQDT] = [] + for nu_ref in nu_ref_list: + nuis = model.calc_channel_nuis(nu_ref) + coefficients = calc_nullvector(model.calc_m_matrix(nu_ref)) + if coefficients is None: + logger.warning("Failed to calculate nullvector for nu_ref=%s, skipping this state.", nu_ref) + continue + coefficients = np.array( + [coeff * (nui ** (3 / 2)) / np.cos(np.pi * nui) for coeff, nui in zip(coefficients, nuis, strict=True)] + ) + coefficients /= np.linalg.norm(coefficients) + arg_max = np.argmax(np.abs(coefficients)) + coefficients *= np.sign(coefficients[arg_max]) + + sqdt_states: list[RydbergStateSQDTAlkalineFJ | RydbergStateSQDTDummy] = [] + for nui, angular_ket in zip(nuis, model.outer_channels, strict=True): + rydberg_class = RydbergStateSQDTDummy if is_dummy_ket(angular_ket) else RydbergStateSQDTAlkalineFJ + sqdt_states.append(rydberg_class.from_angular_ket(model.species, angular_ket, nu=float(nui))) + + states.append(RydbergStateMQDT(coefficients, sqdt_states, nu_ref=nu_ref)) + + return states diff --git a/src/rydstate/utils/__init__.py b/src/rydstate/utils/__init__.py new file mode 100644 index 0000000..e69de29 diff --git a/src/rydstate/utils/linalg.py b/src/rydstate/utils/linalg.py new file mode 100644 index 0000000..cd05b81 --- /dev/null +++ b/src/rydstate/utils/linalg.py @@ -0,0 +1,116 @@ +from __future__ import annotations + +import logging +import math +from typing import TYPE_CHECKING, Any, Literal + +import numpy as np +import scipy +from scipy.optimize import brentq + +if TYPE_CHECKING: + from collections.abc import Callable + from typing import TypeAlias + + import numpy.typing as npt + + NDArray: TypeAlias = npt.NDArray[Any] + +logger = logging.getLogger(__name__) + + +def find_roots( + func: Callable[[float], float], + x_min: float, + x_max: float, + n_grid: int | None = None, + atol: float = 1e-8, +) -> list[float]: + """Find all roots of func in [x_min, x_max]. + + Uses a dense uniform grid to detect sign changes, then refines each bracket + with Brent's method. Validates roots via |func(root)| < atol, which naturally + filters out false sign changes caused by poles (discontinuities of tan). + + Args: + func: 1D scalar function to find roots of. + x_min: Left endpoint of search interval. + x_max: Right endpoint of search interval (must be > x_min). + n_grid: Number of grid points. Defaults to 100 * ceil(x_max - x_min), minimum 1000. + atol: Absolute tolerance for root validation. + + Returns: + Sorted list of x values where func(x) ≈ 0. + + """ + if x_min >= x_max: + return [] + + if n_grid is None: + n_grid = max(50 * math.ceil(x_max - x_min), 500) + + xs = np.linspace(x_min, x_max, n_grid) + fs = np.array([func(x) for x in xs]) + # TODO instead of grid search we also could try to implement a "recursive" kind of method, + # but probably if this is fast enough we dont need it anyway ... + + # Vectorized bracket detection: sign changes between finite adjacent pairs + finite = np.isfinite(fs) + both_finite = finite[:-1] & finite[1:] + sign_change = np.sign(fs[:-1]) * np.sign(fs[1:]) < 0 + bracket_indices = np.where(both_finite & sign_change)[0] + + roots: list[float] = [] + for i in bracket_indices: + try: + root = brentq(func, xs[i], xs[i + 1], xtol=1e-13, rtol=1e-13) + except ValueError: + logger.warning("Brent's method failed to find root in [%f, %f], skipping.", xs[i], xs[i + 1]) + continue + + val = func(root) + if abs(val) > atol: + if abs(val) > 1e5: + pass + # logger.debug("Root not close to zero (probably due to a pole): x=%f f(x)=%e. Skipping.", root, val) # noqa: ERA001, E501 + else: + logger.warning("Root not close to zero: x=%f f(x)=%e. Skipping.", root, val) + continue + + roots.append(root) + + return roots + + +def calc_nullvector( + matrix: NDArray, + *, + method: Literal["numpy_svd", "scipy_nullspace", "scipy_nullspace_gesvd"] = "scipy_nullspace", + atol: float = 1e-8, +) -> NDArray | None: + """Calculate the nullspace vector of a matrix. + + We use scipy.linalg.null_space. + If the nullspace has more than one vector, we raise an error since this should not happen for the MQDT M-matrix. + """ + if matrix.shape == (1, 1): + if abs(matrix[0, 0]) > atol: + raise RuntimeError("Matrix is 1x1 but not close to zero (value=%e), this should not happen.", matrix[0, 0]) + return np.array([1.0]) + + if method == "numpy_svd": + _u, s, vt = np.linalg.svd(matrix) + null_mask = s <= atol + nullspace = vt.T[:, null_mask] + elif method == "scipy_nullspace": + nullspace = scipy.linalg.null_space(matrix, rcond=atol) + elif method == "scipy_nullspace_gesvd": + nullspace = scipy.linalg.null_space(matrix, rcond=atol, lapack_driver="gesvd") + + if nullspace.shape[1] == 0: + logger.error("Nullspace is empty, no solution found.") + return None + if nullspace.shape[1] > 1: + logger.error("Nullspace has more than one vector (shape=%s), returning first vector.", nullspace.shape) + + return np.array(nullspace[:, 0]) From 6bfbc8542e787278e1387eed91282436803e3843 Mon Sep 17 00:00:00 2001 From: johannes-moegerle Date: Tue, 10 Mar 2026 13:46:35 +0000 Subject: [PATCH 10/18] various small fixes --- src/rydstate/angular/angular_ket_dummy.py | 18 ++++++++-- src/rydstate/angular/core_ket_base.py | 6 ++-- src/rydstate/basis/basis_mqdt.py | 5 +-- src/rydstate/rydberg/rydberg_dummy.py | 7 ++-- src/rydstate/rydberg/rydberg_mqdt.py | 36 ++++++++++++------- src/rydstate/species/mqdt/fmodel.py | 13 +++---- .../species/mqdt/species_object_mqdt.py | 4 +-- src/rydstate/utils/linalg.py | 4 +-- 8 files changed, 63 insertions(+), 30 deletions(-) diff --git a/src/rydstate/angular/angular_ket_dummy.py b/src/rydstate/angular/angular_ket_dummy.py index 07775cd..b472979 100644 --- a/src/rydstate/angular/angular_ket_dummy.py +++ b/src/rydstate/angular/angular_ket_dummy.py @@ -35,7 +35,16 @@ def __init__( f_tot: float, m: float | NotSet = NotSet, ) -> None: - """Initialize the Spin ket.""" + """Initialize a dummy angular ket for MQDT perturber channels with unknown quantum numbers. + + Args: + name: Identifier for this dummy channel. By convention the name follows + ``"nl