From 9eccd57d851c10863cc4fcadd9cfb62f4cdcd192 Mon Sep 17 00:00:00 2001 From: Valentin Boettcher Date: Fri, 18 Mar 2022 15:34:49 +0100 Subject: [PATCH] add one qubit model --- one_qubit_model.py | 344 ++++++++++++++++++++++++++++++++++++++++ test_one_qubit_model.py | 38 +++++ 2 files changed, 382 insertions(+) create mode 100644 one_qubit_model.py create mode 100644 test_one_qubit_model.py diff --git a/one_qubit_model.py b/one_qubit_model.py new file mode 100644 index 0000000..ee5a841 --- /dev/null +++ b/one_qubit_model.py @@ -0,0 +1,344 @@ +r""" + Operators for a general model of one qubit coupled to a single bath. + The energy scale is the characteristic energy of the qubit :math:`ω = + 1`. + + The total hamiltonian has the form + + .. math:: + + H=\frac{1}{2}σ_z + + \sqrt{δ} ∑_λ (L^† g_λ b^i_λ + L g_λ^\ast b^{i,†}_λ) + H_B. + + The BCF is normalized so that the integral over its imaginary part is + :math:`-1`. The bath coupling strength is divided by :math:`\langle + L L^†\rangle 2` with respect to the inital state to normalize the + interaction energy to about the order of :math:`ω=1`. + """ + + +from dataclasses import dataclass, field +from numpy.typing import NDArray +from typing import Any, Optional, SupportsFloat +import hops.util.bcf +import hops.util.bcf_fits +import hops.core.hierarchy_parameters as params +import numpy as np +import qutip as qt +from hops.util.abstract_truncation_scheme import TruncationScheme_Simplex +from hops.util.truncation_schemes import ( + TruncationScheme_Power_multi, + BathMemory, +) +import stocproc as sp + +from beartype import beartype +from utility import StocProcTolerances +from model_base import Model +import scipy.special + + +@beartype +@dataclass(eq=False) +class QubitModel(Model): + """ + A class to dynamically calculate all the one qubit model parameters and + generate the HOPS configuration. + + All attributes can be changed after initialization. + """ + + δ: SupportsFloat = 0.1 + """The bath coupling factor.""" + + ω_c: SupportsFloat = 2 + """The cutoff frequency :math:`ω_c`.""" + + s: SupportsFloat = 1 + """The BCF s parameter.""" + + L: qt.Qobj = field(default_factory=lambda: qt.Qobj([[0.0, 1.0], [1.0, 0.0]])) + """ + The :math:`L` coupling operator with shape ``(2, 2)``. + """ + + T: SupportsFloat = 0 + """The temperature of the bath.""" + + ########################################################################### + # HOPS Parameters # + ########################################################################### + + description: str = "" + """A free-form description of the model instance.""" + + bcf_terms: int = 5 + """How many bcf terms to use in the expansions of the BCF.""" + + ψ_0: qt.Qobj = qt.basis([2], [1]) + """The initial state.""" + + t_max: SupportsFloat = 10 + """The maximum simulation time.""" + + resolution: SupportsFloat = 0.1 + """The time resolution of the simulation result.""" + + k_fac: SupportsFloat = 1.7 + """The k_fac parameters for the truncation scheme. + + See + :any:`hops.util.truncation_schemes.TruncationScheme_Power_multi`. + """ + + k_max: int = 10 + """The kmax parameter for the truncation scheme. + + See + :any:`hops.util.abstract_truncation_scheme.TruncationScheme_Simplex` + """ + + influence_tolerance: SupportsFloat = 1e-2 + """The ``influecne_tolerance`` parameter for the truncation + scheme. + + See :any:`hops.util.truncation_schemes.BathMemory`. + """ + + truncation_scheme: str = "bath_memory" + """The truncation scheme to use.""" + + solver_args: dict[str, Any] = field(default_factory=dict) + """Extra arguments for :any:`scipy.integrate.solve_ivp`.""" + + driving_process_tolerance: StocProcTolerances = field( + default_factory=lambda: StocProcTolerances() + ) + """ + The integration and interpolation tolerance for the driving + processes. + """ + + thermal_process_tolerance: StocProcTolerances = field( + default_factory=lambda: StocProcTolerances() + ) + """ + The integration and interpolation tolerance for the thermal noise + processes. + """ + + @property + def system(self) -> qt.Qobj: + """The system hamiltonian.""" + + return 1 / 2 * qt.sigmaz() # type: ignore + + @property + def bcf_norm(self) -> float: + """The normalization factor for the BCF. + + It is not used when generating the stochastic process due to + numerical reasons. It is being incorporated into the + :any:`bcf_scale`. + """ + + return ( + np.pi + * float(self.s) + / ( + scipy.special.gamma(float(self.s) + 1) + * float(self.ω_c) ** float(self.s) + ) + ) + + @property + def bcf_scale(self) -> float: + """ + The BCF scaling factor of the BCF. + """ + + eval = qt.expect(self.L * self.L.dag(), self.ψ_0) + assert isinstance(eval, float) + + return float(self.δ) / (2 * eval.real) * self.bcf_norm + + @property + def bcf(self) -> hops.util.bcf.OhmicBCF_zeroTemp: + """ + The normalized zero temperature BCF. + """ + + return hops.util.bcf.OhmicBCF_zeroTemp( + s=self.s, eta=1, w_c=self.ω_c, normed=False + ) + + @property + def spectral_density(self) -> hops.util.bcf.OhmicSD_zeroTemp: + """ + The normalized zero temperature spectral density. + """ + + return hops.util.bcf.OhmicSD_zeroTemp( + s=float(self.s), + w_c=float(self.ω_c), + eta=1, + normed=False, + ) + + @property + def thermal_correlations( + self, + ) -> Optional[hops.util.bcf.Ohmic_StochasticPotentialCorrelations]: + """ + The normalized thermal noise corellation function. + """ + + if self.T == 0: + return None + + return hops.util.bcf.Ohmic_StochasticPotentialCorrelations( + s=self.s, + eta=1, + w_c=self.ω_c, + normed=False, + beta=1 / float(self.T), + ) + + @property + def thermal_spectral_density( + self, + ) -> Optional[hops.util.bcf.Ohmic_StochasticPotentialDensity]: + """ + The normalized thermal noise spectral density. + """ + + if self.T == 0: + return None + + return hops.util.bcf.Ohmic_StochasticPotentialDensity( + s=self.s, + eta=1, + w_c=self.ω_c, + normed=False, + beta=1.0 / float(self.T), + ) + + def bcf_coefficients( + self, n: Optional[int] = None + ) -> tuple[NDArray[np.complex128], NDArray[np.complex128]]: + """ + The normalizedzero temperature BCF fit coefficients + :math:`G_i,W_i` with ``n`` terms. + """ + + n = n or self.bcf_terms + return self.bcf.exponential_coefficients(n) + + @staticmethod + def basis(n: int = 1) -> qt.Qobj: + """ + A state with of the qubit in the state state ``n`` where ``1`` + means down and ``0`` means up. + """ + + return qt.basis([2], [n]) + + @property + def driving_process(self) -> sp.StocProc: + """The driving stochastic process of the ``i``th bath.""" + + return sp.StocProc_FFT( + spectral_density=self.spectral_density, + alpha=self.bcf, + t_max=self.t_max, + intgr_tol=self.driving_process_tolerance.integration, + intpl_tol=self.driving_process_tolerance.interpolation, + negative_frequencies=False, + ) + + @property + def thermal_process(self) -> Optional[sp.StocProc]: + """The thermal noise stochastic process of the ``i``th bath.""" + + if self.T == 0: + return None + + return sp.StocProc_TanhSinh( + spectral_density=self.thermal_spectral_density, + alpha=self.thermal_correlations, + t_max=self.t_max, + intgr_tol=self.thermal_process_tolerance.integration, + intpl_tol=self.thermal_process_tolerance.interpolation, + negative_frequencies=False, + ) + + ########################################################################### + # Utility # + ########################################################################### + + @property + def hops_config(self): + """ + The hops :any:`hops.core.hierarchy_params.HIParams` parameter object + for this system. + """ + + g, w = self.bcf_coefficients(self.bcf_terms) + + system = params.SysP( + H_sys=self.system.full(), + L=[self.L.full()], + g=[g], + w=[w], + bcf_scale=[self.bcf_scale], + T=[self.T], + description=self.description, + psi0=self.ψ_0.full().flatten(), + ) + + trunc_scheme = TruncationScheme_Power_multi.from_g_w( + g=system.g, + w=system.w, + p=[1, 1], + q=[0.5, 0.5], + kfac=[float(self.k_fac)], + ) + + if self.truncation_scheme == "bath_memory": + trunc_scheme = BathMemory.from_system( + system, + nonlinear=True, + influence_tolerance=float(self.influence_tolerance), + ) + + if self.truncation_scheme == "simplex": + trunc_scheme = TruncationScheme_Simplex(self.k_max) + + hierarchy = params.HiP( + seed=0, + nonlinear=True, + terminator=False, + result_type=params.ResultType.ZEROTH_AND_FIRST_ORDER, + accum_only=False, + rand_skip=None, + truncation_scheme=trunc_scheme, + save_therm_rng_seed=True, + auto_normalize=True, + ) + + default_solver_args = dict(rtol=1e-8, atol=1e-8) + default_solver_args.update(self.solver_args) + + integration = params.IntP( + t_max=float(self.t_max), + t_steps=int(float(self.t_max) / float(self.resolution)) + 1, + **default_solver_args, + ) + + return params.HIParams( + SysP=system, + IntP=integration, + HiP=hierarchy, + Eta=[self.driving_process], + EtaTherm=[self.thermal_process], + ) diff --git a/test_one_qubit_model.py b/test_one_qubit_model.py new file mode 100644 index 0000000..fb4509a --- /dev/null +++ b/test_one_qubit_model.py @@ -0,0 +1,38 @@ +import pytest +import random +from one_qubit_model import * +from qutip import * +from utility import assert_serializable +import scipy.integrate + + +def test_sd_bcf_norm(): + random.seed(0) + + for _ in range(4): + model = QubitModel( + ω_c=random.uniform(0.5, 4), + s=random.choice([1.0, 0.1, 0.5, 0.3]), + T=random.uniform(0.1, 4), + δ=1, + ) + assert_serializable(model) + + bcf = model.bcf + assert np.abs( + scipy.integrate.quad(lambda t: bcf(t).imag, 0, np.inf)[0] * model.bcf_norm + ) == pytest.approx(1) + assert np.sum(model.bcf_coefficients()[0]) == pytest.approx(bcf(0), 1e-2) + + tsd = model.thermal_spectral_density + assert tsd is not None + assert pytest.approx(tsd(1)) == model.spectral_density(1) * 1 / np.expm1( + 1 / float(model.T) + ) + + tcorr = model.thermal_correlations + assert tcorr is not None + + # tests if the normalization of tcorr is correct by + # calculating the fourier transform + model.thermal_process