pyqed package

Subpackages

Submodules

pyqed.Floquet module

pyqed.FranckCondon module

pyqed.backup module

pyqed.beam module

pyqed.cavity module

pyqed.common module

pyqed.coordinates module

pyqed.correlation module

pyqed.fft module

pyqed.fft.dft(x, f, k)

Discrete Fourier transfrom at specified momentum

pyqed.fft.dft2(x, y, f, kx, ky)

Discrete Fourier transfrom at specified momentum

pyqed.fft.fft(a, x=None, axis=-1, **kwargs)

customized 1D Fourier transform of function f along a chosen axis based on numpy

\[g(\omega) = \int dt f(t) * e^{- i * \omega * t}\]
Parameters:

  • f (ndarray) – input data

  • x (TYPE, optional) – the grid points. If None, set to arange(N). The default is None.

  • axis (int, optional) – Axis over which to compute the FFT. If not given, the last axis is used.

  • **kwargs (TYPE) – DESCRIPTION.

Returns:

  • g (ndarray) – the fourier transform of f

  • freq (1darray) – frequencies where f are evaluated

pyqed.fft.fft2(f, dx=1, dy=1)

customized FFT for 2D function input:

f: 2d array,

input array

Returns:

1d array

frequencies

g: 2d array

fourier transform of f

Return type:

freq

pyqed.fft.ifft(a, x=None, axis=-1)

customized fourier transform of function f g = int dt f(t) * exp(i * freq * t) :returns: frequencies where f are evaluated

g: the fourier transform of f

Return type:

freq

pyqed.group module

pyqed.gwp module

pyqed.liouville module

pyqed.mol module

Created on Tue Jun 30 21:16:53 2020

@author: Bing Gu

Basic module for many-level systems

class pyqed.mol.LVC(E, modes)

Bases: Mol

linear vibronic coupling model in Fock space

APES(x)
add_coupling(coupling)

add additional coupling terms to the Hamiltonian such as Stark and Zeeman effects

Parameters:

coupling (list, [[a, b], strength]) – describe the coupling, a, b labels the electronic states

Returns:

updated H.

Return type:

ndarray

buildH()

Calculate the vibronic Hamiltonian.

Parameters:

nums (list of integers) – size for the Fock space of each mode

Returns:

Hamiltonian

Return type:

2d array

buildop(i, f=None, isherm=True)

construct electronic operator

ket{f}ra{i}

if isherm:

return ket{f}ra{i} + ket{i}ra{f}

Parameters:

  • i (int) – initial state.

  • f (int, optional) – final state. Default None. If None, set f = i.

  • isherm (bool, optional) – indicator of whether the returned matrix is Hermitian or not Default: True

Returns:

DESCRIPTION.

Return type:

2darray

calc_edip()
coordinate(n)

build coordinate operators in the full space

Parameters:

n (int) – mode id

Returns:

DESCRIPTION.

Return type:

TYPE

dpes(q)
groundstate()

return the ground state

Returns:

DESCRIPTION.

Return type:

TYPE

plot_PES(x, y)

plot the 3D adiabatic potential energy surfaces

Parameters:

  • x (TYPE) – DESCRIPTION.

  • y (TYPE) – DESCRIPTION.

Returns:

DESCRIPTION.

Return type:

TYPE

promote(A, which='el')
rdm_el(psi)

Compute the electronic reduced density matrix.

Parameters:

psi (TYPE) – DESCRIPTION.

Returns:

DESCRIPTION.

Return type:

TYPE

vertical(n=1)

generate the initial state created by vertical excitation

Parameters:

n (int, optional) – initially excited state. The default is 1.

Returns:

psi – DESCRIPTION.

Return type:

TYPE

wavepacket_dynamics(method='RK4')
class pyqed.mol.Mode(omega: float, couplings: list, truncate: int = 2)

Bases: object

couplings: list
omega: float
truncate: int = 2
class pyqed.mol.Mol(H, edip=None, lowering=None, edip_rms=None, gamma=None)

Bases: object

DQC()
ETPA()
Floquet(omegad, E0, nt)
PE(pump, probe, t2=0.0, **kwargs)

alias for photon_echo

PE2(omega1, omega2, t3=0.0, **kwargs)

2D photon echo signal at -k1+k2+k3 Transforming t1 and t2 to frequency domain.

Parameters:

  • pump (TYPE) – DESCRIPTION.

  • probe (TYPE) – DESCRIPTION.

  • t2 (TYPE, optional) – DESCRIPTION. The default is 0.0.

  • **kwargs (TYPE) – DESCRIPTION.

Returns:

DESCRIPTION.

Return type:

TYPE

TPA()
absorption(omegas, method='sos', **kwargs)

Linear absorption of the model.

Parameters:

  • omegas (TYPE) – DESCRIPTION.

  • method (TYPE, optional) – DESCRIPTION. The default is ‘SOS’.

  • normalize (TYPE, optional) – DESCRIPTION. The default is True.

Returns:

DESCRIPTION.

Return type:

TYPE

cars(shift, omega1, t2=0.0, plt_signal=False, fname=None)
deom(bath, coupling, coupling_dipole=None, pulse_system_func=None, pulse_coupling_func=None, mode=None)

hierarchical equations of motion

driven_dynamics(psi0, pulse, dt=0.001, Nt=1, obs_ops=None, nout=1, t0=0.0)

wavepacket dynamics in the presence of laser pulses

Parameters:

  • pulse (TYPE) – DESCRIPTION.

  • dt (TYPE, optional) – DESCRIPTION. The default is 0.001.

  • Nt (TYPE, optional) – DESCRIPTION. The default is 1.

  • obs_ops (TYPE, optional) – DESCRIPTION. The default is None.

  • nout (TYPE, optional) – DESCRIPTION. The default is 1.

  • t0 (float) – initial time

Return type:

None.

property edip
property edip_rms
eigenenergies()
eigenstates(k=None)
Parameters:

k (integer) – number of eigenstates to compute, < dim

Returns:

  • eigvals (vector)

  • eigvecs (2d array)

eigvals()
energy(psi)
evolve(psi0=None, dt=0.001, Nt=10, e_ops=None, nout=1, t0=0.0, pulse=None)

quantum dynamics under time-independent Hamiltonian

Parameters:

  • pulse (TYPE) – DESCRIPTION.

  • dt (TYPE, optional) – DESCRIPTION. The default is 0.001.

  • Nt (TYPE, optional) – DESCRIPTION. The default is 1.

  • obs_ops (TYPE, optional) – DESCRIPTION. The default is None.

  • nout (TYPE, optional) – DESCRIPTION. The default is 1.

  • t0 (float) – initial time

Return type:

None.

getH()
get_dip()
get_dm()
get_edip()
get_nonhermH()
get_nonhermitianH()
get_p_from_r()
\[p_{ji} = -i m (\omega_i - \omega_j) x_{ji} = i \omega_{ij} \mu_ji\]
Return type:

None.

ground_state(method='trivial')
groundstate(method='trivial')
multi(nmol=2)
photon_echo(pump, probe, t2=0.0, **kwargs)

2D photon echo signal at -k1+k2+k3

Parameters:

  • pump (TYPE) – DESCRIPTION.

  • probe (TYPE) – DESCRIPTION.

  • t2 (TYPE, optional) – DESCRIPTION. The default is 0.0.

  • **kwargs (TYPE) – DESCRIPTION.

Returns:

DESCRIPTION.

Return type:

TYPE

quantum_dynamics(psi0, dt=0.001, Nt=1, obs_ops=None, nout=1, t0=0.0)

quantum dynamics under time-independent hamiltonian

Parameters:

  • pulse (TYPE) – DESCRIPTION.

  • dt (TYPE, optional) – DESCRIPTION. The default is 0.001.

  • Nt (TYPE, optional) – DESCRIPTION. The default is 1.

  • obs_ops (TYPE, optional) – DESCRIPTION. The default is None.

  • nout (TYPE, optional) – DESCRIPTION. The default is 1.

  • t0 (float) – initial time

Return type:

None.

run(psi0=None, dt=0.01, e_ops=None, nt=1, nout=1, t0=0.0, edip=None, pulse=None)

quantum dynamics under time-independent Hamiltonian

Parameters:

  • pulse (TYPE) – DESCRIPTION.

  • dt (float, optional) – time interval. The default is 0.001.

  • Nt (int, optional) – int. The default is 1.

  • obs_ops (TYPE, optional) – DESCRIPTION. The default is None.

  • nout (TYPE, optional) – DESCRIPTION. The default is 1.

  • t0 (float) – initial time

Return type:

None.

set_decay(gamma)

Set the decay rate for all excited states.

Parameters:

gamma (TYPE) – DESCRIPTION.

Return type:

None.

set_decay_for_all(gamma)

Set the decay rate for all excited states.

Parameters:

gamma (TYPE) – DESCRIPTION.

Return type:

None.

set_dephasing(gamma)

set the pure dephasing rate for all coherences

Parameters:

gamma (TYPE) – DESCRIPTION.

Return type:

None.

set_dipole(dip)
set_edip(edip, pol=None)
set_lifetime(tau)
set_mdip(mdip)
tcl2()

second-order time-convolutionless quantum master equation

class pyqed.mol.Result(description=None, psi0=None, rho0=None, dt=None, Nt=None, times=None, t0=0, nout=1)

Bases: object

analyze()
dump(fname)

save results to disk

Parameters:

fname (TYPE) – DESCRIPTION.

Return type:

None.

expect()
save(fname)
class pyqed.mol.SESolver(H=None)

Bases: object

correlation_2op_1t()
correlation_3op_1t(psi0, oplist, dt, Nt)

<AB(t)C>

Parameters:

  • psi0 –

  • oplist –

  • dt –

  • Nt –

correlation_3op_2t(psi0, oplist, dt, Nt, Ntau)

<A(t)B(t+tau)C(t)> :param oplist: [a, b, c] :type oplist: list of arrays :param psi0: initial state :type psi0: array :param dt: :param nt: number of time steps for t :type nt: integer :param ntau: time steps for tau :type ntau: integer

correlation_4op_1t(psi0, oplist, dt=0.005, Nt=1)

<AB(t)C(t)D>

Parameters:

  • psi0 –

  • oplist –

  • dt –

  • Nt –

correlation_4op_2t(psi0, oplist, dt=0.005, Nt=1, Ntau=1)
Parameters:

  • psi0 (vector) – initial state

  • oplist (list of arrays) –

propagator(dt, Nt)
run(psi0=None, dt=0.01, Nt=1, e_ops=None, nout=1, t0=0.0, edip=None, pulse=None, use_sparse=True)

quantum dynamics under time-independent and time-dependent Hamiltonian

Parameters:

  • psi0 (1d array) – initial state

  • pulse (callable/list of callable) – external laser pulses

  • dt (TYPE, optional) – DESCRIPTION. The default is 0.001.

  • Nt (TYPE, optional) – DESCRIPTION. The default is 1.

  • obs_ops (TYPE, optional) – DESCRIPTION. The default is None.

  • nout (TYPE, optional) – DESCRIPTION. The default is 1.

  • t0 (float) – initial time. The default is 0.

Return type:

None.

pyqed.mol.direct_product(a, n)

Single-molecule operator in the direct product space of N molecules

\[A = a \otimes I \cdots I + I \otimes A \otimes I \cdots I + \cdots\]
Parameters:

  • a (TYPE) – DESCRIPTION.

  • n (TYPE) – DESCRIPTION.

Returns:

atot – DESCRIPTION.

Return type:

TYPE

pyqed.mol.driven_dynamics(H, psi0, dt=0.01, Nt=1, e_ops=None, nout=1, t0=0.0, return_result=True, use_sparse=True)

Laser-driven dynamics in the presence of laser pulses

Parameters:

  • ham (2d array) – Hamiltonian of the molecule

  • dip (TYPE) – transition dipole moment

  • psi0 (1d array) – initial wavefunction

  • pulse (TYPE) – laser pulse

  • dt (TYPE) – time step.

  • Nt (TYPE) – timesteps.

  • e_ops (list) – observable operators to compute

  • nout (int) – Store data every nout steps

  • sparse (bool) – indicater whether to use sparse matrices or not (is there a gain?)

Return type:

None.

pyqed.mol.high_frequency_drive()
pyqed.mol.load_result(fname)

read result obj saved with pickle

pyqed.mol.mls(dim=3)
pyqed.mol.polar(x, y)
pyqed.mol.read_input(fname_e, fname_edip, g_included=True)

Read input data from quantum chemistry output.

Parameters:

  • fname_e (str) – filename for the energy levels.

  • fname_edip (list) – filenames for the electric dipole moment.

  • g_included (TYPE, optional) – DESCRIPTION. The default is True.

Returns:

mol – DESCRIPTION.

Return type:

TYPE

pyqed.mol.test_sesolver()

pyqed.noise module

pyqed.optics module

Created on Tue Mar 26 17:26:02 2019

@author: binggu

class pyqed.optics.Analyser(E, t)

Bases: object

FROG(w=None, use_fft=False)
plot_spectrogram()
propagator()
spectrogram()
class pyqed.optics.Biphoton(omegap, bw, Te, p=None, q=None, phase_matching='sinc')

Bases: object

bandwidth(which='signal')

Compute the bandwidth of the signal/idler mode

Parameters:

which (TYPE, optional) – DESCRIPTION. The default is ‘signal’.

Return type:

None.

detect()

two-photon detection amplitude in a temporal grid defined by the spectral grid.

Returns:

  • t1 (1d array)

  • t2 (1d array)

  • d (detection amplitude in the temporal grid (t1, t2))

detect_is()
detect_si()
g2()
get_jsa()
Returns:

jsa – joint spectral amplitude

Return type:

array

get_jta()

Compute the joint temporal amplitude J(ts, ti) over a temporal meshgrid.

Returns:

  • ts (1d array) – signal time grid

  • ti (1d array) – idler temporal grid

  • jta (2d array) – joint temporal amplitude

jta(ts, ti)
plt_jsa(xlabel=None, ylabel=None, fname=None)
pump(bandwidth)

pump pulse envelope :param bandwidth:

rdm(which='signal')
set_grid(p, q)
class pyqed.optics.ChirpedPulse(omegac=0.11024710050125681, tau=206.70686668237954, tc=0, amplitude=0.001, cep=0.0, beta=0)

Bases: Pulse

efield(t)
Parameters:

t (TYPE) – DESCRIPTION.

Return type:

electric field at time t.

spectrum(omega)

Fourier transform of the Gaussian pulse

class pyqed.optics.GaussianPulse(omegac=0.11024710050125681, tau=206.70686668237954, tc=0, delay=0.0, amplitude=0.001, cep=0.0, beta=0, polarization=None)

Bases: object

efield(t)
Parameters:

t (TYPE) – DESCRIPTION.

Return type:

electric field at time t.

envelop(t)
field(t)

electric field

plt_efield()
spectrogram(efield)
spectrum(omega)

Fourier transform of the Gaussian pulse

class pyqed.optics.Mode(freq, nmax, decay=None, polarization=None)

Bases: object

annihilate()
create()
getH(ZPE=False)
get_annihilate()
get_create()
get_nonhermH()
get_nonhermitianH()

non-Hermitian Hamiltonian for the cavity mode

Params:

kappa: decay constant

Returns:

DESCRIPTION.

Return type:

TYPE

get_num()

number operator

ham(ZPE=False)
number_operator()

number operator input:

N: integer

number of states

setQ(Q)
vacuum(sparse=True)

get initial density matrix for cavity vacuum state

vacuum_dm()

get initial density matrix for cavity vacuum state

class pyqed.optics.Pulse(omegac=0.11024710050125681, tau=206.70686668237954, tc=0, delay=0.0, amplitude=0.001, intensity=None, cep=0.0, beta=0, polarization=None)

Bases: object

E(t)
efield(t, return_complex=False)
Parameters:

t (TYPE) – time.

Returns:

complex-valued electric field at time t. Take the real part for the electric field.

Return type:

complex

envelop(t)
field(t)

electric field

plt_efield()
spectrogram(efield)
spectrum(omega)

Fourier transform of the Gaussian pulse

pyqed.optics.field_to_intensity(E)
pyqed.optics.fwhm_to_std(fwhm)
pyqed.optics.hom(p, q, f, tau)

HOM coincidence probability

Parameters:

  • p –

  • q –

  • f –

  • tau –

  • method (str) –

    “brute”: directly integrating the JSA over the frequency grid “schmidt”: compute the signal using the Schmidt modes of the

    entangled light

  • nmodes –

Returns:

prob – coincidence probability

Return type:

1d array

pyqed.optics.hom_schmidt(p, q, f, method='rdm', nmodes=5)

HOM signal with Schmidt modes

Parameters:

  • p –

  • q –

  • f –

  • nmodes –

pyqed.optics.intensity_to_field(I)

transform intensity to electric field

E = sqrt(

rac{2I}{c epsilon_0})

ITYPE

intensity, in units of W/cm^2

None.

electric field in atomic units

pyqed.optics.jta(t2, t1, omegap, sigmap, Te)

Analytical form for the joint temporal amplitude for SPDC type-II two-photon state.

Note that two single-photon electric field prefactors are neglected.

Parameters:

  • t2 (TYPE) – DESCRIPTION.

  • t1 (TYPE) – DESCRIPTION.

Return type:

None.

pyqed.optics.rdm(f, dx=1, dy=1, which='x')

Compute the reduced density matrix by tracing out the other dof for a 2D wavefunction

Parameters:

  • f (2D array) – 2D wavefunction

  • dx (float, optional) – DESCRIPTION. The default is 1.

  • dy (float, optional) – DESCRIPTION. The default is 1.

  • which (str) – indicator which rdm is required. Default is ‘x’.

Returns:

rho1 – Reduced density matrix

Return type:

TYPE

pyqed.optics.schmidt_decompose(f, dp, dq, nmodes=5, method='rdm')

kernel method f: 2D array,

input function to be decomposed

nmodes: int

number of modes to be kept

method: str

rdm or svd

pyqed.optics.std_to_fwhm(tau)

Transformation between standard deviation to FWHM for a Gaussian pulse

Parameters:

tau (float) – std

Returns:

FWHM.

Return type:

float

pyqed.oqs module

pyqed.phys module

class pyqed.phys.HarmonicOscillator(omega, mass=1, x0=0)

Bases: object

basic class for harmonic oscillator

eigenstate(x, n=0)
potential(x)
class pyqed.phys.HeisenbergModel(Jx, Jy, Jz, h, nsites)

Bases: object

DMRG()
NARG()
build()
run()
class pyqed.phys.Morse(D, a, re, mass=1)

Bases: object

basic class for Morse oscillator

eigenstate(x, n=0)
eigval(n)

Given an (integer) quantum number v, return the vibrational energy.

Parameters:

n (TYPE) – DESCRIPTION.

Return type:

None.

potential(x)
class pyqed.phys.TFIM(J, g, nsites)

Bases: object

build()
renormalize()
run(nroots=1)
pyqed.phys.anticomm(A, B)
pyqed.phys.anticommutator(A, B)
pyqed.phys.basis(N, j)
Parameters:

  • N (int) – Size of Hilbert space for a multi-level system.

  • j (int) – The j-th basis function.

Returns:

j-th basis function for the Hilbert space.

Return type:

1d complex array

pyqed.phys.basis_transform(A, v)

transformation rule: A_{ab} = <a|i><i|A|j><j|b> = Anew = v^dag A v input:

A: matrix of operator A in old basis v: basis transformation matrix

output:

Anew: matrix A in the new basis

pyqed.phys.boson(omega, n, ZPE=False)
pyqed.phys.cartesian(*args)

compute Cartesian product of args

pyqed.phys.cartesian2polar(x, y)

transform Cartesian coordinates to polar

\[x = r \cos( heta) y = r \sin( heta)\]
Parameters:

  • r (TYPE) – DESCRIPTION.

  • theta (TYPE) – DESCRIPTION.

Returns:

  • x (TYPE) – DESCRIPTION.

  • y (TYPE) – DESCRIPTION.

pyqed.phys.cartesian_product(arrays)

A fast cartesion product function

pyqed.phys.coh_op(j, i, d)

return a matrix representing the coherence :math:`|j

angle langle i|`

jTYPE

DESCRIPTION.

iTYPE

DESCRIPTION.

dTYPE

DESCRIPTION.

None.

pyqed.phys.comm(A, B)
pyqed.phys.commutator(A, B)
pyqed.phys.coth(x)
pyqed.phys.dag(a)
pyqed.phys.dagger(a)
pyqed.phys.destroy(N)

Annihilation operator for bosons.

Parameters:

N (int) – Size of Hilbert space.

Return type:

2d array complex

pyqed.phys.discretize(a=0, b=1, l=4, endpoints=True)

Create a uniform math with size with level l in the range [a, b]

mesh size is \((b-a)/2^l\)

Parameters:

  • a (TYPE, optional) – DESCRIPTION. The default is 0.

  • b (TYPE, optional) – DESCRIPTION. The default is 1.

  • l (TYPE, optional) – DESCRIPTION. The default is 4.

  • endpoints (TYPE, optional) – DESCRIPTION. The default is False.

Returns:

DESCRIPTION.

Return type:

TYPE

pyqed.phys.driven_dissipative_dynamics(ham, dip, rho0, pulse, dt=0.001, Nt=1, obs_ops=None, nout=1)

Laser-driven dynamics in the presence of laser pulses

Parameters:

  • ham (2d array) – Hamiltonian of the molecule

  • dip (TYPE) – transition dipole moment

  • rho0 (2d array complex) – initial density matrix

  • pulse (TYPE) – laser pulse

  • dt (float) – DESCRIPTION.

  • Nt (TYPE) – DESCRIPTION.

  • obs_ops (list) – observable operators to compute

  • nout (int) – Store data every nout steps

Return type:

None.

pyqed.phys.driven_dynamics(ham, dip, psi0, pulse, dt=0.001, Nt=1, obs_ops=None, nout=1, t0=0.0)

Laser-driven dynamics in the presence of laser pulses

Parameters:

  • ham (2d array) – Hamiltonian of the molecule

  • dip (TYPE) – transition dipole moment

  • psi0 (1d array) – initial wavefunction

  • pulse (TYPE) – laser pulse

  • dt (TYPE) – time step.

  • Nt (TYPE) – timesteps.

  • obs_ops (list) – observable operators to compute

  • nout (int) – Store data every nout steps

Return type:

None.

pyqed.phys.eig_asymm(h)

Diagonalize a real, asymmetrix matrix and return sorted results.

Return the eigenvalues and eigenvectors (column matrix) sorted from lowest to highest eigenvalue.

pyqed.phys.eigh(a, k=None, **args)

find the eigenvalues and eigenstates of matrix a

Parameters:

  • a (TYPE) – DESCRIPTION.

  • k (TYPE, optional) – DESCRIPTION. The default is None.

Returns:

DESCRIPTION.

Return type:

TYPE

pyqed.phys.expect(psi, op, variance=False)

compute the expectation value and its variances

\[\langle \psi | O | \psi\]

angle

psiTYPE

DESCRIPTION.

opTYPE

DESCRIPTION.

varianceTYPE, optional

DESCRIPTION. The default is False.

TYPE

DESCRIPTION.

pyqed.phys.expm(A, t, method='EOM')
exponentiate a matrix at t

U(t) = e^{A t}

Parameters:

  • A (TYPE) – DESCRIPTION.

  • t (float or list) – times

  • method (TYPE, optional) –

    DESCRIPTION. The default is ‘EOM’.

    EOM: equation of motion approach.

    d/dt U(t) = A U(t) This can be generalized for time-dependent Hamiltonians A(t)

    diagonalization: diagonalize A

    for Hermitian matrices only, this is prefered

Returns:

Ulist – DESCRIPTION.

Return type:

TYPE

pyqed.phys.fermi(E, Ef=0.0, T=0.0001)

Fermi-Dirac distribution function INPUT:

E : Energy Ef : Fermi energy T : temperture (in units of energy, i.e., kT)

OUTPUT:

f(E): Fermi-Dirac distribution function at energy E

pyqed.phys.fftfreq(times)

get the spectral range corresponding to the temporal grid (a.u.)

Parameters:

times (TYPE) – DESCRIPTION.

Returns:

DESCRIPTION.

Return type:

TYPE

pyqed.phys.gaussian(x, sigma=1.0)

normalized Gaussian distribution

Returns:

DESCRIPTION.

Return type:

TYPE

pyqed.phys.get_index(array, value)

get the index of element in array closest to value

pyqed.phys.gwp(x, a=None, x0=0.0, p0=0.0, ndim=1)

complex Gaussian wavepacket

\[g(x; x_0, p_0) = Det(A)^{1/4}/\pi^{n/4} e^{-1/2 (x-x_0) A (x-x_0) + i p_0(x-x_0)}\]
Parameters:

  • x (TYPE) – DESCRIPTION.

  • sigma (TYPE, optional) – (co)variance matrix.

  • x0 (TYPE, optional) – DESCRIPTION. The default is 0..

  • p0 (TYPE, optional) – DESCRIPTION. The default is 0..

Returns:

psi – DESCRIPTION.

Return type:

TYPE

pyqed.phys.gwp2(x, y, sigma=array([[1., 0.], [0., 1.]]), xc=[0, 0], kc=[0, 0])

generate a 2D Gaussian wavepacket in grid :param x0: float, mean value of gaussian wavepacket along x :param y0: float, mean value of gaussian wavepacket along y :param sigma: float array, covariance matrix with 2X2 dimension :param kx0: float, initial momentum along x :param ky0: float, initial momentum along y :return: float 2darray, the gaussian distribution in 2D grid

pyqed.phys.gwp_k(k, sigma, x0, k0)

analytical fourier transform of gauss_x(x), above

pyqed.phys.ham_ho(freq, n, ZPE=False)

Hamiltonian for harmonic oscilator

input:

freq: fundemental frequency in units of Energy n : size of matrix ZPE: boolean, if ZPE is included in the Hamiltonian

output:

h: hamiltonian of the harmonic oscilator

pyqed.phys.heaviside(x)
pyqed.phys.hilbert_dist(A, B)

Returns the Hilbert-Schmidt distance between two density matrices A & B.

Parameters:

  • A (qobj) – Density matrix or state vector.

  • B (qobj) – Density matrix or state vector with same dimensions as A.

Returns:

dist – Hilbert-Schmidt distance between density matrices.

Return type:

float

Notes

See V. Vedral and M. B. Plenio, Phys. Rev. A 57, 1619 (1998).

pyqed.phys.integrate(f, a, b, **args)

Compute a definite integral.

Parameters:

  • f (TYPE) – DESCRIPTION.

  • a (TYPE) – DESCRIPTION.

  • b (TYPE) – DESCRIPTION.

  • **args (TYPE) – DESCRIPTION.

Returns:

DESCRIPTION.

Return type:

TYPE

pyqed.phys.interval(x)
pyqed.phys.is_positive_def(A)
pyqed.phys.isdiag(M)

Check if a matrix is diagonal.

Parameters:

M (TYPE) – DESCRIPTION.

Returns:

DESCRIPTION.

Return type:

TYPE

pyqed.phys.isherm(a)
pyqed.phys.isunitary(m)
pyqed.phys.jacobi_anger(n, z=1)

Jacobi-Anger expansion .. math:

e^{iz \cos(     heta)} = \sum_{n=-\infty}^\infty i^n J_n(z) e^{in       heta}
Parameters:

  • n (TYPE) – DESCRIPTION.

  • z (TYPE, optional) – DESCRIPTION. The default is 1.

Returns:

DESCRIPTION.

Return type:

TYPE

pyqed.phys.jump(f, i, dim=2, isherm=True)
pyqed.phys.ket2dm(psi)

convert a ket into a density matrix

Parameters:

psi (TYPE) – DESCRIPTION.

Returns:

DESCRIPTION.

Return type:

TYPE

pyqed.phys.ldo(b, A)

linear differential operator Ab

Parameters:

  • b (TYPE) – DESCRIPTION.

  • A (TYPE) – DESCRIPTION.

Returns:

DESCRIPTION.

Return type:

TYPE

pyqed.phys.lindbladian(l, rho)

lindblad superoperator: l rho l^dag - 1/2 * {l^dag l, rho} l is the operator corresponding to the disired physical process e.g. l = a, for the cavity decay and l = sm for polarization decay

pyqed.phys.liouvillian(rho, H, c_ops)

lindblad quantum master eqution

pyqed.phys.logarithmic_discretize(n, base=2)

log discretization of (0, 1) .. math:

[\Lambda^{-(n+1)}, \Lambda^{-n}], n = 0, 1, ...
Parameters:

  • n (TYPE) – DESCRIPTION.

  • base (TYPE, optional) – discretization parameter. The default is 2.

Returns:

\(\Lambda^{-n}, n = 0, 1 ...\) in descending order.

Return type:

TYPE

pyqed.phys.lorentzian(x, width=1.0)
Parameters:

  • x (TYPE) – DESCRIPTION.

  • x0 (float) – center of the Lorentzian

  • width (float) – Half-wdith half-maximum

Return type:

None.

pyqed.phys.lowering(dims=2)
pyqed.phys.meshgrid(*args)

fix the indexing of the Numpy meshgrid

Parameters:

*args (TYPE) – DESCRIPTION.

Returns:

DESCRIPTION.

Return type:

TYPE

pyqed.phys.morse(r, D, a, re)
Morse potential

D * (1. - e^{-a * (r - r_e)})**2

Parameters:

  • r (float/1darray) – DESCRIPTION.

  • D (float) – well depth

  • a (TYPE) – ‘width’ of the potential.

  • re (TYPE) – equilibrium bond distance

Returns:

DESCRIPTION.

Return type:

TYPE

pyqed.phys.multi_spin(onsite, nsites)

construct the hamiltonian for a multi-spin system params:

onsite: array, transition energy for each spin nsites: number of spins

Returns:

  • ham (ndarray) – Hamiltonian

  • lower (ndnarry) – lowering operator

pyqed.phys.multiboson(omega, nmodes, J=0, truncate=2)

construct the hamiltonian for a multi-boson system

Parameters:

  • omegas (1D array) – resonance frequenies of the boson modes

  • nmodes (integer) – number of boson modes

  • J (float) – hopping constant

  • truncation (TYPE, optional) – DESCRIPTION. The default is 2.

Returns:

  • ham (TYPE) – DESCRIPTION.

  • lower (TYPE) – DESCRIPTION.

pyqed.phys.multimode(omegas, nmodes, J=0, truncate=2)

construct the direct tensor-product Hamiltonian for a multi-mode system

Parameters:

  • omegas (1D array) – resonance frequenies of the boson modes

  • nmodes (integer) – number of boson modes

  • J (float) – nearest-neighbour hopping constant

  • truncate (list) – size of Fock space for each mode

Returns:

  • ham (TYPE) – DESCRIPTION.

  • xs (list) – position operators in the composite space for each mode

pyqed.phys.multispin(onsite, hopping, nsites)
pyqed.phys.nlargest(a, n=1, with_index=False)

finds the largest n elements from a Python iterable

Parameters:

  • a (TYPE) – DESCRIPTION.

  • n (int) – DESCRIPTION.

  • with_index (bool) – if index is needed.

Returns:

DESCRIPTION.

Return type:

TYPE

pyqed.phys.norm(psi, dx=1)

normalization of the wavefunction

N = int dx psi^*(x) psi(x)

Parameters:

psi (1d array, complex) – DESCRIPTION.

Return type:

float, L2 norm

pyqed.phys.norm2(f, dx=1, dy=1)

L2 norm of the 2D array f

Parameters:

  • f (TYPE) – DESCRIPTION.

  • dx (TYPE) – DESCRIPTION.

  • dy (TYPE) – DESCRIPTION.

Returns:

DESCRIPTION.

Return type:

TYPE

pyqed.phys.obs(psi, a)
Parameters:

  • psi (1d array) – wavefunction.

  • a (2d array) – operator a.

Returns:

Expectation of operator a.

Return type:

complex

pyqed.phys.obs_dm(rho, d)

observables for operator d

pyqed.phys.overlap(bra, ket)
pyqed.phys.pauli()
pyqed.phys.pdf_normal(x, mu=0, sigma=1.0)
pyqed.phys.polar2cartesian(r, theta)

transform polar coordinates to Cartesian

\[x = r \cos( heta) y = r \sin( heta)\]
Parameters:

  • r (TYPE) – DESCRIPTION.

  • theta (TYPE) – DESCRIPTION.

Returns:

  • x (TYPE) – DESCRIPTION.

  • y (TYPE) – DESCRIPTION.

pyqed.phys.polarization_vector(pol='x')

unit length polarization vector

Parameters:

d (TYPE, optional) – DESCRIPTION. The default is ‘x’.

Returns:

DESCRIPTION.

Return type:

TYPE

pyqed.phys.project(P, a)

reduce the representation of operators to a subspace defined by the projection operator P

Parameters:

  • P (TYPE) – DESCRIPTION.

  • a (TYPE) – DESCRIPTION.

Return type:

None.

pyqed.phys.propagator(H, dt, nt)

compute the propagator for time-dependent and time-independent H U(t) = e^{-i H t}

Parameters:

  • H (ndarray or list of ndarray or callable) – Hamiltonian. if H is ndarray, H is time-independent. Otherwise H is time-dependent.

  • t (float or list) – times

pyqed.phys.propagator_H_const(H, dt, nt, method='diag')

compute the propagator for time-dependent and time-independent H U(t) = e^{-i H t}

Parameters:

  • H (ndarray or list of ndarray or callable) – Hamiltonian. if H is ndarray, H is time-independent. Otherwise H is time-dependent.

  • t (float or list) – times

pyqed.phys.ptrace(rho, dims, which='B')

partial trace of subsystems in a density matrix defined in a composite space

Parameters:

  • rho (ndarray) – DESCRIPTION.

  • which (TYPE, optional) – DESCRIPTION. The default is ‘B’.

Returns:

rhoA – DESCRIPTION.

Return type:

TYPE

pyqed.phys.quadrature(n)

Quadrature operator of a photon mode X = (a + a+)/sqrt{2}

Parameters:

n (TYPE) – DESCRIPTION.

Returns:

DESCRIPTION.

Return type:

TYPE

pyqed.phys.quantum_dynamics(ham, psi0, dt=0.001, Nt=1, obs_ops=None, nout=1, t0=0.0, output='obs.dat')

Laser-driven dynamics in the presence of laser pulses

Parameters:

  • ham (2d array) – Hamiltonian of the molecule

  • psi0 (1d array) – initial wavefunction

  • dt (float) – time step.

  • Nt (int) – timesteps.

  • obs_ops (list) – observable operators to compute

  • nout (int) – Store data every nout steps

Return type:

None.

pyqed.phys.raising(dims=2)

raising operator for spin-1/2 :param dims: Hilbert space dimension :type dims: integer

Returns:

sp – raising operator

Return type:

2x2 array

pyqed.phys.rect(x)
pyqed.phys.resolvent(omega, Ulist, dt)

compute the resolvent 1/(omega - H) from the Fourier transform of the propagator omega: float

frequency to evaluate the resolvent

Ulist: list of matrices

propagators

dt: time-step used in the computation of U

pyqed.phys.rgwp(x, x0=0.0, sigma=1.0)

real Gaussian wavepacket

Parameters:

  • x (TYPE) – DESCRIPTION.

  • x0 (float) – central position

  • sigma (TYPE) – DESCRIPTION.

Return type:

None.

pyqed.phys.rk4(rho, fun, dt, *args)

Runge-Kutta method

pyqed.phys.rotate(angle)
pyqed.phys.sigmax()
pyqed.phys.sigmay()
pyqed.phys.sigmaz()
pyqed.phys.sinc(x)
\[\]

sinc(x) =

rac{sin(x)}{x}

xTYPE

DESCRIPTION.

TYPE

DESCRIPTION.

pyqed.phys.sort(eigvals, eigvecs)

sort eigenvalues and eigenvectors from low to high

Parameters:

  • eigvals (TYPE) – DESCRIPTION.

  • eigvecs (TYPE) – DESCRIPTION.

Returns:

  • eigvals (TYPE) – DESCRIPTION.

  • eigvecs (TYPE) – DESCRIPTION.

pyqed.phys.spin_ops(m)

Spin operators represented by Sz eigenstates

Parameters:

m (int) – spin multiplicity \(m = 2S + 1\)

Raises:

NotImplementedError – DESCRIPTION.

Returns:

DESCRIPTION.

Return type:

TYPE

pyqed.phys.stepsize(x)
pyqed.phys.tdse(wf, h)
pyqed.phys.tensor(*args, **kwargs)

Calculates the tensor product of input operators.

Build from QuTip.

Parameters:

args (array_like) – list or array of quantum objects for tensor product.

pyqed.phys.tensor_power(a, n: int)

kron(a, kron(a, …))

pyqed.phys.thermal_dm(n, u)

return the thermal density matrix for a boson n: integer

dimension of the Fock space

u: float

reduced temperature, omega/k_B T

pyqed.phys.tracedist(A, B)

Calculates the trace distance between two density matrices.. See: Nielsen & Chuang, “Quantum Computation and Quantum Information”

Parameters ———-!= A : 2D array (N,N)

Density matrix or state vector.

B2D array (N,N)

Density matrix or state vector with same dimensions as A.

tolfloat

Tolerance used by sparse eigensolver, if used. (0=Machine precision)

sparse{False, True}

Use sparse eigensolver.

Returns:

tracedist – Trace distance between A and B.

Return type:

float

Examples

>>> x=fock_dm(5,3)
>>> y=coherent_dm(5,1)
>>> tracedist(x,y)
0.9705143161472971
pyqed.phys.transform(A, v)

transformation rule: A_{ab} = <a|i><i|A|j><j|b> = Anew = v^dag A v input:

A: matrix of operator A in old basis v: basis transformation matrix v[old_index, new_index]

output:

Anew: matrix A in the new basis

pyqed.qnm module

pyqed.quadrature module

pyqed.spin_boson module

pyqed.style module

pyqed.style.color_code(x, y, z, fig, ax, cbar=False)
pyqed.style.curve(x, y, **kwargs)

simple 1D curve plot

Parameters:

  • x (TYPE) – DESCRIPTION.

  • y (TYPE) – DESCRIPTION.

Return type:

None.

pyqed.style.export(x, y, z, fname='output.dat', fmt='gnuplot')

export 3D data to gnuplot format

Parameters:

  • x (TYPE) – DESCRIPTION.

  • y (TYPE) – DESCRIPTION.

  • z (TYPE) – DESCRIPTION.

  • fname (TYPE, optional) – DESCRIPTION. The default is ‘output.dat’.

  • fmt (str, optional) – The target format. The default is ‘gnuplot’.

Return type:

None.

pyqed.style.imshow(x, y, f, vmin=None, vmax=None, ticks=None, output='output.pdf', xlabel='X', ylabel='Y', diverge=False, cmap='viridis', **kwargs)
Parameters:

  • f (2D array) – array to be plotted.

  • extent (list [xmin, xmax, ymin, ymax]) –

Return type:

Save a fig in the current directory.

pyqed.style.level_scheme(E, ylim=None, fname=None)

plot the energy levels :param E:

pyqed.style.matplot(x, y, f, vmin=None, vmax=None, ticks=None, output='output.pdf', xlabel='X', ylabel='Y', diverge=False, cmap='viridis', **kwargs)
Parameters:

  • f (2D array) – array to be plotted.

  • extent (list [xmin, xmax, ymin, ymax]) –

Returns:

  • Save a fig in the current directory.

  • To be deprecated. Please use imshow.

pyqed.style.plot_surface(x, y, surface)
pyqed.style.plot_surfaces(x, y, surfaces)
pyqed.style.read_result(fname)

read result obj saved with pickle

pyqed.style.scatter(points)
pyqed.style.set_style(fontsize=12)
pyqed.style.subplots(nrows=1, ncols=1, figsize=(4, 3), sharex=True, sharey=True, **kwargs)
pyqed.style.surf(x, y, f, fname='output.png', xlabel='X', ylabel='Y', zlabel='Z', cmap=None, title=None, method='matplotlib')
pyqed.style.test_level_scheme()
pyqed.style.tocolor(x, vmin=0, vmax=1)

map values to colors

Parameters:

x (TYPE) – DESCRIPTION.

Returns:

rbga color.

Return type:

TYPE

pyqed.style.two_scales(x, yl, yr, xlabel=None, ylabels=None, xlim=None, yllim=None, yrlim=None, yticks=None, fname='output.pdf')
pyqed.style.vector_field(f, **kwargs)

3D plot of a vector field.

Parameters:

f (list of x,y,z components) – DESCRIPTION.

Return type:

None.

pyqed.superoperator module

pyqed.susceptibility module

pyqed.units module

class pyqed.units.AtomicUnits

Bases: object

pyqed.wigner module

Created on Sat Jul 24 22:23:15 2021

@author: bing

@Source: https://www.frank-zalkow.de/en/the-wigner-ville-distribution-with-python.html

pyqed.wigner.nextpow2(p)
pyqed.wigner.spectrogram(x, d=1)

Wigner transform of an input signal with FFT.

W(w, t) = int d au x(t + au/2) x^*(t - au/2) e^{i w tau}

Parameters:

  • x (1d array) – The time-domain signal.

  • d (TYPE, optional) – DESCRIPTION. The default is 1.

Returns:

  • TYPE – spectrogram in (f, t)

  • freqs (1d array) – sample frequencies.

pyqed.wigner.wigner(x, d=1)

Wigner transform of an input signal with FFT.

For a state vector

\[W(t, w) = \int d au x(t + au/2) x^*(t - au/2) e^{i w au}\]

For density matrices

\[W(t, w) = \int d au\]

ho(t + tau/2, t - tau/2) e^{i w au}

xTYPE

a state vector/density matrix.

TYPE

DESCRIPTION.

freqs

ranges of the frequency/momentum space

pyqed.wigner.wvd(audioFile, t=None, N=None, trace=0, make_analytic=True)

pyqed.wpd module

Created on Mon Jan 4 23:44:15 2021

Exact nonadiabatic wavepacket dynamics solver for vibronic models with N vibrational modes (N = 1 ,2)

For linear coordinates, use SPO method For curvilinear coordinates, use RK4 method

@author: Bing Gu

pyqed.wpd.KEO(psi, kx, ky, G)

compute kinetic energy operator K * psi

Parameters:

  • psi (TYPE) – DESCRIPTION.

  • dt (TYPE) – DESCRIPTION.

Return type:

None.

pyqed.wpd.PEO(psi, v)

V |psi> :param dt: float

time step

Parameters:

  • v_2d – float array the two electronic states potential operator in grid basis

  • psi_grid – list the two-electronic-states vibrational state in grid basis

Returns:

psi_grid(update): list the two-electronic-states vibrational state in grid basis after being half time step forward

class pyqed.wpd.ResultSPO2(**args)

Bases: Result

dump(fname)

save results to disk

Parameters:

fname (TYPE) – DESCRIPTION.

Return type:

None.

get_population(fname=None, plot=False)
plot_population(p)
plot_wavepacket(psilist, state_id=None, **kwargs)
position(plot=False, fname=None)
pyqed.wpd.S0(x, y)
pyqed.wpd.S1(x, y)
class pyqed.wpd.SPO(x, mass=1, nstates=1)

Bases: object

Time-dependent Schrodinger Equation for wavepackets on a single PES with 1D nuclear coordinate

build(dt)
k_evolve(psi_x)

one time step for exp(-i * K * dt)

run(psi0, dt, nt=1, t0=0, nout=1)

Time-dependent Schrodinger Equation for wavepackets on a single PES.

Parameters:

  • psi0 (1d array, complex) – initial wavepacket

  • t0 (float) – initial time

  • dt (float) – the small time interval over which to integrate

  • nt (float, optional) – the number of intervals to compute. The total change in time at the end of this method will be dt * Nsteps. default is N = 1

set_grid(xmin=-1, xmax=1, npts=32)
set_potential(potential)
x_evolve(psi)

one time step for exp(-i * V * dt)

Parameters:

psi (TYPE) – DESCRIPTION.

Returns:

DESCRIPTION.

Return type:

TYPE

x_evolve_half(psi)

one time step for exp(-i * V * dt)

Parameters:

psi (TYPE) – DESCRIPTION.

Returns:

DESCRIPTION.

Return type:

TYPE

class pyqed.wpd.SPO2(x, y, mass=None, nstates=2, coords='linear', G=None, abc=False)

Bases: object

second-order split-operator method for nonadiabatic wavepacket dynamics in the diabatic representation with two-dimensional nuclear coordinate

For time-independent Hamiltonian

e^{-i H Delta t} = e^{- i V Delta t/2} e^{-i K Delta t} e^{-iVDelta t/2}

For time-dependent H,

TBI

build(dt, inertia=None)

Setup the propagators appearing in the split-operator method.

For the kinetic energy operator with Jacobi coordinates

K =

rac{p_r^2}{2mu} + rac{1}{I(r)} p_ heta^2

Since the two KEOs for each dof do not commute, it has to be factorized as

e^{-i K delta t} = e{-i K_1 delta t} e^{- i K_2 delta t}

where $p_ heta = -i pa_ heta$ is the momentum operator.

dtTYPE

DESCRIPTION.

intertia: func

moment of inertia, only used for Jocabi coordinates.

None.

current_density(psi, state_id=0)

Compute the velocity field of the vibrational flow

Parameters:

  • psi (TYPE) – DESCRIPTION.

  • state_id (TYPE, optional) – DESCRIPTION. The default is 0.

Return type:

None.

plot_surface(style='2D')
plt_wp(psilist, **kwargs)
population(psi, representation='diabatic', plot=False)

return the electronic populations

Parameters:

psi (TYPE) – DESCRIPTION.

Return type:

2darray, nt, nstates

rdm_el(psi)

Compute the reduced electronic density matrix

Parameters:

psi (TYPE) – vibronic wavefunction.

Returns:

rho – DESCRIPTION.

Return type:

TYPE

run(psi0, e_ops=[], dt=0.01, nt=1, t0=0.0, nout=1, return_states=True)
setG(G)
set_DPES(surfaces, diabatic_couplings, eta=None)

set the potential energy operatpr from the diabatic PES and vibronic couplings

To be deprecated.

Parameters:

  • surfaces (TYPE) – DESCRIPTION.

  • diabatic_couplings (TYPE) – DESCRIPTION.

  • abc (boolean, optional) – indicator of whether using absorbing boundary condition. This is often used for dissociation. The default is False.

Returns:

v – DESCRIPTION.

Return type:

TYPE

set_dpes(v)
set_grid(x, y)
set_masses(mass)
class pyqed.wpd.SPO2NH(x, y, *args, **kwargs)

Bases: SPO2

build(dt)

Setup the propagators appearing in the split-operator method.

For the kinetic energy operator with Jacobi coordinates

K =

rac{p_r^2}{2mu} + rac{1}{I(r)} p_ heta^2

Since the two KEOs for each dof do not commute, it has to be factorized as

e^{-i K delta t} = e{-i K_1 delta t} e^{- i K_2 delta t}

where $p_ heta = -i pa_ heta$ is the momentum operator.

dtTYPE

DESCRIPTION.

intertia: func

moment of inertia, only used for Jocabi coordinates.

None.

position(psilist, plot=False)
run(psi0, e_ops=[], dt=0.01, nt=1, t0=0.0, nout=1, return_states=True)
class pyqed.wpd.SPO3(x, y, z, masses, nstates=2, coords='linear', G=None, abc=False)

Bases: object

second-order split-operator method for nonadiabatic wavepacket dynamics in the diabatic representation with three-dimensional nuclear space

For time-independent Hamiltonian

\[e^{-i H \Delta t} = e^{- i V \Delta t/2} e^{-i K \Delta t} e^{-iV\Delta t/2}\]
For time-dependent H,

# TODO

build(dt, inertia=None)

Setup the propagators appearing in the split-operator method.

For the kinetic energy operator with Jacobi coordinates

K =

rac{p_r^2}{2mu} + rac{1}{I(r)} p_ heta^2

Since the two KEOs for each dof do not commute, it has to be factorized as

e^{-i K delta t} = e{-i K_1 delta t} e^{- i K_2 delta t}

where $p_ heta = -i pa_ heta$ is the momentum operator.

dtTYPE

DESCRIPTION.

intertia: func

moment of inertia, only used for Jocabi coordinates.

None.

plot_surface(style='2D')
plt_wp(psilist, **kwargs)
population(psi)

return the electronic populations

Parameters:

psi (TYPE) – DESCRIPTION.

Return type:

None.

run(psi0, e_ops=[], dt=0.01, nt=1, t0=0.0, nout=1, return_states=True)
setG(G)
set_DPES(surfaces, diabatic_couplings, eta=None)

set the diabatic PES and vibronic couplings

Parameters:

  • surfaces (TYPE) – DESCRIPTION.

  • diabatic_couplings (TYPE) – DESCRIPTION.

  • abc (boolean, optional) – indicator of whether using absorbing boundary condition. This is often used for dissociation. The default is False.

Returns:

v – DESCRIPTION.

Return type:

TYPE

set_abc()
set_grid(x, y, z)
set_masses(masses)
class pyqed.wpd.Solver

Bases: object

set_obs_ops(obs_ops)
pyqed.wpd.adiabatic_2d(x, y, psi0, v, dt, Nt=0, coords='linear', mass=None, G=None)

propagate the adiabatic dynamics at a single surface

Parameters:

  • dt – time step

  • v – 2d array potential matrices in 2D

  • psi – list the initial state

mass: list of 2 elements

reduced mass

Nt: int
the number of the time steps, Nt=0 indicates that no propagation has been done,

only the initial state and the initial purity would be the output

G: 4D array nx, ny, ndim, ndim

G-matrix

Returns:

psi_end: list the final state

G: 2d array

G matrix only used for curvilinear coordinates

pyqed.wpd.diabatic(x, y)

PESs in diabatic representation :param x_range_half: float, the displacement of potential from the origin

in x

Parameters:

  • y_range_half – float, the displacement of potential from the origin in y

  • couple_strength – the coupling strength between these two potentials

  • couple_type – int, the nonadiabatic coupling type. here, we used: 0) no coupling 1) constant coupling 2) linear coupling

Returns:

v: float 2d array, matrix elements of the DPES and couplings

pyqed.wpd.diabatic_coupling(x, y)
pyqed.wpd.divergence(f, h)

div(F) = dFx/dx + dFy/dy + … g = np.gradient(Fx,dx, axis=1)+ np.gradient(Fy,dy, axis=0) #2D g = np.gradient(Fx,dx, axis=2)+ np.gradient(Fy,dy, axis=1) +np.gradient(Fz,dz,axis=0) #3D

pyqed.wpd.dpsi(psi, kx, ky, ndim=2)

Momentum operator operates on the wavefunction

Parameters:

  • psi (2D complex array) – DESCRIPTION.

  • ndim (int, default 2) – coordinates dimension

Returns:

kpsi – DESCRIPTION.

Return type:

(nx, ny, ndim)

pyqed.wpd.dxpsi(psi)

Momentum operator operates on the wavefunction

Parameters:

psi (2D complex array) – DESCRIPTION.

Returns:

kpsi – DESCRIPTION.

Return type:

(nx, ny, ndim)

pyqed.wpd.dypsi(psi)

Momentum operator operates on the wavefunction

Parameters:

psi (2D complex array) – DESCRIPTION.

Returns:

kpsi – DESCRIPTION.

Return type:

(nx, ny, ndim)

pyqed.wpd.gauss_k(k, a, x0, k0)

analytical fourier transform of gauss_x(x), above

pyqed.wpd.gauss_x_2d(sigma, x0, y0, kx0, ky0)

generate the gaussian distribution in 2D grid :param x0: float, mean value of gaussian wavepacket along x :param y0: float, mean value of gaussian wavepacket along y :param sigma: float array, covariance matrix with 2X2 dimension :param kx0: float, initial momentum along x :param ky0: float, initial momentum along y :return: gauss_2d: float array, the gaussian distribution in 2D grid

pyqed.wpd.hpsi(psi, kx, ky, v, G)
pyqed.wpd.k_evolve_2d(dt, masses, kx, ky, psi)

propagate the state in grid basis a time step forward with H = K :param dt: float, time step :param kx: float, momentum corresponding to x :param ky: float, momentum corresponding to y :param psi_grid: list, the two-electronic-states vibrational states in

grid basis

Returns:

psi_grid(update): list, the two-electronic-states vibrational states in grid basis

pyqed.wpd.plot_wavepacket(x, y, psilist, **kwargs)
pyqed.wpd.potential_2d(x_range_half, y_range_half, couple_strength, couple_type)

generate two symmetric harmonic potentials wrt the origin point in 2D :param x_range_half: float, the displacement of potential from the origin

in x

Parameters:

  • y_range_half – float, the displacement of potential from the origin in y

  • couple_strength – the coupling strength between these two potentials

  • couple_type – int, the nonadiabatic coupling type. here, we used: 0) no coupling 1) constant coupling 2) linear coupling

Returns:

v_2d: float list, a list containing for matrices: v_2d[0]: the first potential matrix v_2d[1]: the potential coupling matrix

between the first and second

v_2d[2]: the potential coupling matrix

between the second and first

v_2d[3]: the second potential matrix

pyqed.wpd.square_barrier(x, width, height)
pyqed.wpd.theta(x)
theta function :

returns 0 if x<=0, and 1 if x>0

pyqed.wpd.x_evolve_2d(dt, psi, v)

propagate the state in grid basis half time step forward with H = V :param dt: float

time step

Parameters:

  • v_2d – float array the two electronic states potential operator in grid basis

  • psi_grid – list the two-electronic-states vibrational state in grid basis

Returns:

psi_grid(update): list the two-electronic-states vibrational state in grid basis after being half time step forward

Module contents