pytdscf.basis package

Submodules

pytdscf.basis.abc module

Abstract DVR primitive class module

class pytdscf.basis.abc.DVRPrimitivesMixin(ngrid: int)[source]

Bases: ABC

Abstract DVR primitive function class

abstract diagnalize_pos_rep_matrix() → None[source]

Numerical diagonalization of pos_rep_matrix.

If analytical diagonalization is known, implemented in inheritance class

dvr_func(n: int, q: float) → float[source]

DVR function

\[\chi_\alpha=\sum_{j=0}^{N-1}\varphi_j(q)U_{j\alpha} \quad (\alpha=0,\ldots, N-1)\]

In other words,

\[|\chi_\alpha\rangle =U^\dagger |\varphi_j\rangle\]

abstract fbr_func(n: int, q: float) → float[source]: fbrinal Primitive function \(\varphi_n(q)\), such as HO, Sine, etc

get_1st_derivative_matrix_dvr() → ndarray[source]: \(\langle\chi_\alpha|\frac{d}{dq}|\chi_\beta\rangle\)

abstract get_1st_derivative_matrix_fbr() → ndarray[source]

Numerical integral of \(\langle\varphi_j|\frac{d}{dq}|\varphi_k\rangle\)

If analytical integral is known, implemented in inheritance class

get_2nd_derivative_matrix_dvr() → ndarray[source]: \(\langle\chi_\alpha|\frac{d^2}{dq^2}|\chi_\beta\rangle\)

abstract get_2nd_derivative_matrix_fbr() → ndarray[source]

Numerical integral of \(\langle\varphi_j|\frac{d^2}{dq^2}|\varphi_k\rangle\)

If analytical integral is known, implemented in inheritance class

get_grids() → List[float][source]: grids \(x_\alpha\) correspond to eigenvalue of pos_rep_matrix

abstract get_pos_rep_matrix() → ndarray[source]

Numerical integral of \(\langle\varphi_j|\hat{q}|\varphi_k\rangle\)

If analytical integral is known, implemented in inheritance class

get_sqrt_weights(k: int = 0) → List[float][source]: \(\sqrt{w_\alpha}=U_{k\alpha}^{\ast}/\varphi_k(x_\alpha)\)

get_unitary() → ndarray[source]

Get Unitary Matrix which diagonalize pos_rep_matrix

Returns:: u[alpha,j] = \((U_{j\alpha})^\dagger\) = \((U^\dagger)_{\alpha j}\)
Return type:: np.ndarray

where,

\[\sum_{j,k} U_{j\alpha}\langle\varphi_j|\hat{q}|\varphi_k\rangle U_{k\beta}^\dagger = x_\alpha \delta_{\alpha\beta}\]

plot_dvr(n: int = None, q: int = None) → None[source]: Plot DVR functions \(\{\chi_n(q)\}\)

plot_fbr(n: int = None, q: float = None) → None[source]: Plot FBR \(\{\varphi_n(q)\}\)

pytdscf.basis.boson module

Boson basis module

class pytdscf.basis.boson.Boson(nstate: int)[source]

Bases: object

Boson basis class

Parameters:: nstate (int) – number of states

get_annihilation_matrix(margin: int = 0) → ndarray[tuple[int, ...], dtype[float64]][source]

Returns the annihilation matrix for the boson basis

\[b |n> = \sqrt{n} |n-1>\]

Example

>>> exciton = Boson(4)
>>> exciton.get_annihilation_matrix()
array([[0., 1., 0., 0.],
       [0., 0., 1.41421356, 0.],
       [0., 0., 0., 1.73205081],
       [0., 0., 0., 0.]])
>>> state = np.array([[1.0],
                      [1.0],
                      [1.0],
                      [0.0]])
>>> np.allclose(exciton.get_annihilation_matrix() @ state, np.array([[0.0], [1.0], [math.sqrt(2)], [math.sqrt(3)]]))
True

get_creation_matrix(margin: int = 0) → ndarray[tuple[int, ...], dtype[float64]][source]

Returns the creation matrix for the boson basis

\[b^\dagger |n> = \sqrt{n+1} |n+1>\]

Example

>>> boson = Boson(3)
>>> boson.get_creation_matrix()
array([[0., 0., 0.],
       [1., 0., 0.],
       [1.41421356, 0., 0.]])
>>> state = np.array([[1.0],
                      [1.0],
                      [0.0]])
>>> np.allclose(boson.get_creation_matrix() @ state, np.array([[0.0], [1.0], [math.sqrt(2)]]))
True

get_number_matrix() → ndarray[tuple[int, ...], dtype[float64]][source]

Returns the number matrix for the boson basis

\[b^\dagger b |n> = n |n>\]

Example

>>> boson = Boson(3)
>>> boson.get_number_matrix()
array([[0., 0., 0.],
       [0., 1., 0.],
       [0., 0., 2.]])
>>> state = np.array([[1.0],
                      [1.0],
                      [0.0]])
>>> np.allclose(boson.get_number_matrix() @ state, np.array([[0.0], [1.0], [0.0]]))
True

get_p2_matrix() → ndarray[tuple[int, ...], dtype[float64]][source]: Returns the p^2 matrix for the boson basis

\[p^2 = -1/2 (b^\dagger b + b b^\dagger - b^\dagger ^2 - b^2)\]

get_p_matrix(margin: int = 0) → ndarray[tuple[int, ...], dtype[complex128]][source]: Returns the p matrix for the boson basis

\[p = 1/\sqrt{2} i (b - b^\dagger) = i/\sqrt{2} (b^\dagger - b)\]

get_q2_matrix() → ndarray[tuple[int, ...], dtype[float64]][source]: Returns the p^2 matrix for the boson basis

\[p^2 = 1/2 (b^\dagger b + b b^\dagger + b^\dagger ^2 + b^2)\]

get_q_matrix(margin: int = 0) → ndarray[tuple[int, ...], dtype[float64]][source]

Returns the q matrix for the boson basis

q = 1/sqrt(2) (b + b^dagger)

names: list[str]

property nprim: int

nstate: int

pytdscf.basis.exciton module

Exciton basis module

class pytdscf.basis.exciton.Exciton(nstate: int, names: list[str] | None = None)[source]

Bases: object

Exciton basis class

Parameters:

nstate (int) – number of states
names (list[str], optional) – names of the states such as [“down”, “up”]. Defaults to None, in which case the names are set to [“S0”, “S1”, …].

get_annihilation_matrix() → ndarray[tuple[int, ...], dtype[float64]][source]

Returns the annihilation matrix for the exciton basis

Example

>>> exciton = Exciton(2)
>>> exciton.get_annihilation_matrix()
array([[0., 1.],
       [0., 0.]])
>>> state = np.array([[1.0],
                      [0.0]])
>>> np.allclose(exciton.get_annihilation_matrix() @ state, np.array([[0.0], [0.0]]))
True
>>> state = np.array([[0.0],
                      [1.0]])
>>> np.allclose(exciton.get_annihilation_matrix() @ state, np.array([[1.0], [0.0]]))

get_creation_matrix() → ndarray[tuple[int, ...], dtype[float64]][source]

Returns the creation matrix for the exciton basis

Example

>>> exciton = Exciton(2)
>>> exciton.get_creation_matrix()
array([[0., 0.],
       [1., 0.]])
>>> state = np.array([[1.0],
                      [0.0]])
>>> np.allclose(exciton.get_creation_matrix() @ state, np.array([[0.0], [1.0]]))
True
>>> state = np.array([[0.0],
                      [1.0]])
>>> np.allclose(exciton.get_creation_matrix() @ state, np.array([[0.0], [0.0]]))
True

names: list[str]

property nprim: int

nstate: int

pytdscf.basis.exponential module

Exponential DVR module

class pytdscf.basis.exponential.Exponential(ngrid: int, length: float, x0: float = 0.0, doAnalytical: bool = True)[source]

Bases: DVRPrimitivesMixin

Exponential DVR is equivalent to FFT. It is suitable for periodic system.

pytdscf.basis.ho module

The Harmonic Oscillator (HO) primitive basis functions and DVR functions.

class pytdscf.basis.ho.HarmonicOscillator(ngrid: int, omega: float, q_eq: float | None = 0.0, units: str | None = 'cm-1', dimnsionless: bool | None = False, doAnalytical: bool | None = True)[source]

Bases: DVRPrimitivesMixin

Harmonic Oscillator DVR functions

See also MCTDH review Phys.Rep. 324, 1 (2000) appendix B https://doi.org/10.1016/S0370-1573(99)00047-2

Normalization factor :

\[A_n = \frac{1}{\sqrt{n! 2^n}} \left(\frac{m\omega}{\pi\hbar}\right)^{\frac{1}{4}}\ \xrightarrow[\rm impl]{} \frac{1}{\sqrt{n! 2^n}} \left(\frac{\omega_i}{\pi}\right)^{\frac{1}{4}}\]

Dimensionless coordinate :

\[\zeta = \sqrt{\frac{m\omega} {\hbar}}(x-a) \xrightarrow[\rm impl]{} \sqrt{\omega_i}q_i\]

Primitive Function :

\[\varphi_n = A_n H_n(\zeta) \exp\left(- \frac{\zeta^2}{2}\right)\quad (n=0,2,\ldots,N-1)\]

ngrid

# of grid

Type:: int

nprim

# of primitive function. same as ngrid

Type:: int

omega

frequency in a.u.

Type:: float

q_eq

eq_point in mass-weighted coordinate

Type:: float

dimensionless

Input is dimensionless coordinate or not.

Type:: bool

doAnalytical

Use analytical integral or diagonalization.

Type:: bool

diagnalize_pos_rep_matrix()[source]: Analytical formulation has not yet derived.

fbr_func(n: int, q)[source]

Primitive functions :math:\varphi_n(q)

Parameters:: n (int) – index (0 <= n < ngrid)
Returns:: \(\varphi_n(q)\)
Return type:: float

get_1st_derivative_matrix_dvr()[source]: \(\langle\chi_\alpha|\frac{d}{dq}|\chi_\beta\rangle\)

get_1st_derivative_matrix_fbr()[source]: Analytical Formulation

\[D_{j k}^{(1)}=\ \left\langle\varphi_{j} | \frac{\mathrm{d}}{\mathrm{d} x} | \varphi_{k}\right\rangle=\ -\sqrt{\frac{m \omega}{2}}\left(\sqrt{j+1} \delta_{j, k-1}-\sqrt{j} \delta_{j, k+1}\right)\]

get_2nd_derivative_matrix_dvr()[source]: \(\langle\chi_\alpha|\frac{d^2}{dq^2}|\chi_\beta\rangle\)

get_2nd_derivative_matrix_fbr()[source]: Analytical Formulation

\[D_{j k}^{(2)}=\ \frac{m\omega}{2}\left(\sqrt{(j-1)j}\delta_{j,k+2}-(2j+1)\delta_{j,k} + \sqrt{(j+2)(j+1)}\delta_{j,k-2}\right)\]

get_ovi_CS_HO(p: float, q: float, type: str = 'DVR') → ndarray[source]

Get overlap integral 1D-array between coherent state <p,q,ω,| and |HO(ω)>

Parameters:

p (float) – momentum of coherent state in mass-weighted a.u.
q (float) – position of coherent state in mass-weighted a.u.
type (str) – Whether “DVR” or “FBR”. Default is “DVR”.

Returns:

overlap integral 1D-array between coherent state <p,q,ω,| and |HO(ω)>

Return type:

np.ndarray

get_pos_rep_matrix()[source]: HO has analytical formulation

\[\left\langle\varphi_{j}|\hat{x}| \varphi_{k}\right\rangle=\ \sqrt{\frac{j+1}{2 m \omega}} \delta_{j, k-1}+x_{\mathrm{eq}}\ \delta_{j k}+\sqrt{\frac{j}{2 m \omega}} \delta_{j, k+1}\]

class pytdscf.basis.ho.PrimBas_HO(origin: float, freq_cm1: float, nprim: int, origin_is_dimless: bool = True)[source]

Bases: object

The Harmonic Oscillator eigenfunction primitive basis.

This class holds information on Gauss Hermite type SPF basis functions. The n-th primitive represents

Normalization factor :

\[A_n = \frac{1}{\sqrt{n! 2^n}} \left(\frac{m\omega}{\pi\hbar}\right)^{\frac{1}{4}}\]

Dimensionless coordinate :

\[\zeta = \sqrt{\frac{m\omega} {\hbar}}(x-a)\]

Primitive Function :

\[\chi_n = A_n H_n(\zeta) \exp\left(- \frac{\zeta^2}{2}\right)\]

Parameters:

origin (float) – The center (equilibrium) dimensionless coordinate \(\sqrt{\frac{m\omega}{\hbar}} a\) of Hermite polynomial.
freq_cm1 (float) – The frequency \(\omega\) of Hermite polynomial. The unit is cm-1.
nprim (int) – The number (max order) of primitive basis on a certain SPF.
origin_is_dimless (bool, optional) – If True, given self.origin is dimensionless coordinate. Defaults to True.

origin

The center (equilibrium) dimensionless coordinate \(\sqrt{\frac{m\omega}{\hbar}} a\) of Hermite polynomial.

Type:: float

freq_cm1

The frequency \(\omega\) of Hermite polynomial. The unit is cm-1.

Type:: float

nprim

The number (max order) of primitive basis on a certain SPF.

Type:: int

freq_au

self.freq_cm1 in a.u.

Type:: float

origin_mwc

Mass weighted coordinate self.origin in a.u.

Type:: float

origin: origin is mass-weghted

todvr()[source]

pytdscf.basis.sin module

Sine DVR module

class pytdscf.basis.sin.Sine(ngrid: int, length: float, x0: float = 0.0, units: str = 'angstrom', doAnalytical: bool = True, include_terminal: bool = True)[source]

Bases: DVRPrimitivesMixin

Sine DVR functions Note that Sine DVR position matrix is not tridiagonal! Index starts from j=1, alpha=1

Terminal points (x_{0} and x_{N+1}) do not belog to the grid.

See also MCTDH review Phys.Rep. 324, 1 (2000) appendix B https://doi.org/10.1016/S0370-1573(99)00047-2

Primitive Function :

\[\begin{split}\varphi_{j}(x)= \begin{cases}\sqrt{2 / L} \sin \left(j \pi\left(x-x_{0}\right) / L\right) \ & \text { for } x_{0} \leq x \leq x_{N+1} \\ 0 & \text { else }\end{cases} \quad (j=1,2,\ldots,N)\end{split}\]

ngrid

Number of grid. Note that which does not contain terminal point.

Type:: int

nprim

Number of primitive function. sama as ngrid.

Type:: int

length

Length in a.u.

Type:: float

x0

start point in a.u.

Type:: float

doAnalytical

Use analytical integral or diagonalization.

Type:: bool

include_terminal

Whether include terminal grid.

Type:: bool

diagnalize_pos_rep_matrix()[source]

Analytical diagonalization has been derived

\[U_{j \alpha}=\sqrt{\frac{2}{N+1}} \sin \left(\frac{j \alpha \pi}{N+1}\right)\]

\[z_{\alpha}=\cos \left(\frac{\alpha \pi}{N+1}\right)\]

This leads to the DVR grid points

\[x_{\alpha}=x_{0}+\frac{L}{\pi} \arccos \left(z_{\alpha}\right)\ =x_{0}+\alpha \frac{L}{N+1}=x_{0}+\alpha \Delta x\]

fbr_func(n: int, x: float)[source]

Primitive functions :math:\varphi_n(x)

Parameters:: n (int) – index (0 <= n < ngrid)
Returns:: \(\varphi_n(x)\)
Return type:: float

get_1st_derivative_matrix_dvr()[source]: \(\langle\chi_\alpha|\frac{d}{dq}|\chi_\beta\rangle\)

get_1st_derivative_matrix_fbr()[source]: Analytical Formulation

\[D_{j k}^{(1)}=\ \left\langle\varphi_{j} |\frac{\mathrm{d}}{\mathrm{d} x} | \varphi_{k}\right\rangle={\rm mod} (j-k, 2) \frac{4}{L} \frac{j k}{j^{2}-k^{2}},\ \quad j \neq k\]

get_2nd_derivative_matrix_dvr()[source]: Analytical forumulation exists.

\[\begin{split}D_{\alpha \beta}^{(2), \mathrm{DVR}}=-\left(\frac{\pi}{\Delta x}\right)^{2}\left\{\begin{array}{l} -\frac{1}{3}+\frac{1}{6(N+1)^{2}}-\frac{1}{2(N+1)^{2} \sin ^{2}\left(\frac{\alpha \pi}{N+1}\right)}, \quad \alpha=\beta \\ \frac{2(-1)^{\alpha-\beta}}{(N+1)^{2}} \frac{\sin \left(\frac{\alpha \pi}{N+1}\right) \sin \left(\frac{\beta \pi}{N+1}\right)}{\left(\cos \left(\frac{\alpha \pi}{N+1}\right)-\cos \left(\frac{\beta \pi}{N+1}\right)\right)^{2}},\quad \alpha \neq \beta \end{array} \right.\end{split}\]

get_2nd_derivative_matrix_fbr()[source]: Analytical Formulation (start from j=1 !)

\[D_{j k}^{(1)}=\ \left\langle\varphi_{j} |\frac{\mathrm{d}}{\mathrm{d} x} | \varphi_{k}\right\rangle={\rm mod} (j-k, 2) \frac{4}{L} \frac{j k}{j^{2}-k^{2}},\ \quad j \neq k\]

get_pos_rep_matrix()[source]

Sine position matrix

\[Q_{j k}=\left\langle\varphi_{j}|\hat{z}| \ \varphi_{k}\right\rangle=\frac{1}{2}\left(\delta_{j, k+1}+\delta_{j, k-1}\right)\]

where transformed variable

\[z=\cos \left(\pi\left(x-x_{0}\right) / L\right)\]

is introduced.

Module contents

Caution: The DVR implementation of this module is deprecated. Use the Discvar library instead.

class pytdscf.basis.Boson(nstate: int)[source]

Bases: object

Boson basis class

Parameters:: nstate (int) – number of states

get_annihilation_matrix(margin: int = 0) → ndarray[tuple[int, ...], dtype[float64]][source]

Returns the annihilation matrix for the boson basis

\[b |n> = \sqrt{n} |n-1>\]

Example

>>> exciton = Boson(4)
>>> exciton.get_annihilation_matrix()
array([[0., 1., 0., 0.],
       [0., 0., 1.41421356, 0.],
       [0., 0., 0., 1.73205081],
       [0., 0., 0., 0.]])
>>> state = np.array([[1.0],
                      [1.0],
                      [1.0],
                      [0.0]])
>>> np.allclose(exciton.get_annihilation_matrix() @ state, np.array([[0.0], [1.0], [math.sqrt(2)], [math.sqrt(3)]]))
True

get_creation_matrix(margin: int = 0) → ndarray[tuple[int, ...], dtype[float64]][source]

Returns the creation matrix for the boson basis

\[b^\dagger |n> = \sqrt{n+1} |n+1>\]

Example

>>> boson = Boson(3)
>>> boson.get_creation_matrix()
array([[0., 0., 0.],
       [1., 0., 0.],
       [1.41421356, 0., 0.]])
>>> state = np.array([[1.0],
                      [1.0],
                      [0.0]])
>>> np.allclose(boson.get_creation_matrix() @ state, np.array([[0.0], [1.0], [math.sqrt(2)]]))
True

get_number_matrix() → ndarray[tuple[int, ...], dtype[float64]][source]

Returns the number matrix for the boson basis

\[b^\dagger b |n> = n |n>\]

Example

>>> boson = Boson(3)
>>> boson.get_number_matrix()
array([[0., 0., 0.],
       [0., 1., 0.],
       [0., 0., 2.]])
>>> state = np.array([[1.0],
                      [1.0],
                      [0.0]])
>>> np.allclose(boson.get_number_matrix() @ state, np.array([[0.0], [1.0], [0.0]]))
True

get_p2_matrix() → ndarray[tuple[int, ...], dtype[float64]][source]: Returns the p^2 matrix for the boson basis

\[p^2 = -1/2 (b^\dagger b + b b^\dagger - b^\dagger ^2 - b^2)\]

get_p_matrix(margin: int = 0) → ndarray[tuple[int, ...], dtype[complex128]][source]: Returns the p matrix for the boson basis

\[p = 1/\sqrt{2} i (b - b^\dagger) = i/\sqrt{2} (b^\dagger - b)\]

get_q2_matrix() → ndarray[tuple[int, ...], dtype[float64]][source]: Returns the p^2 matrix for the boson basis

\[p^2 = 1/2 (b^\dagger b + b b^\dagger + b^\dagger ^2 + b^2)\]

get_q_matrix(margin: int = 0) → ndarray[tuple[int, ...], dtype[float64]][source]

Returns the q matrix for the boson basis

q = 1/sqrt(2) (b + b^dagger)

names: list[str]

property nprim: int

nstate: int

class pytdscf.basis.Exciton(nstate: int, names: list[str] | None = None)[source]

Bases: object

Exciton basis class

Parameters:

nstate (int) – number of states
names (list[str], optional) – names of the states such as [“down”, “up”]. Defaults to None, in which case the names are set to [“S0”, “S1”, …].

get_annihilation_matrix() → ndarray[tuple[int, ...], dtype[float64]][source]

Returns the annihilation matrix for the exciton basis

Example

>>> exciton = Exciton(2)
>>> exciton.get_annihilation_matrix()
array([[0., 1.],
       [0., 0.]])
>>> state = np.array([[1.0],
                      [0.0]])
>>> np.allclose(exciton.get_annihilation_matrix() @ state, np.array([[0.0], [0.0]]))
True
>>> state = np.array([[0.0],
                      [1.0]])
>>> np.allclose(exciton.get_annihilation_matrix() @ state, np.array([[1.0], [0.0]]))

get_creation_matrix() → ndarray[tuple[int, ...], dtype[float64]][source]

Returns the creation matrix for the exciton basis

Example

>>> exciton = Exciton(2)
>>> exciton.get_creation_matrix()
array([[0., 0.],
       [1., 0.]])
>>> state = np.array([[1.0],
                      [0.0]])
>>> np.allclose(exciton.get_creation_matrix() @ state, np.array([[0.0], [1.0]]))
True
>>> state = np.array([[0.0],
                      [1.0]])
>>> np.allclose(exciton.get_creation_matrix() @ state, np.array([[0.0], [0.0]]))
True

names: list[str]

property nprim: int

nstate: int

class pytdscf.basis.Exponential(ngrid: int, length: float, x0: float = 0.0, doAnalytical: bool = True)[source]

Bases: DVRPrimitivesMixin

Exponential DVR is equivalent to FFT. It is suitable for periodic system.