API / Code Reference

Submodules

pytalises.wavefunction module

The Wavefunction class and its attributes.

class pytalises.wavefunction.Wavefunction(psi, number_of_grid_points, spatial_ext, t0=0.0, m=1.054571817e-34, variables={}, normalize_const=None)

Bases: object

Class describing wave function.

Parameters:
  • psi (string list of strings) – Each string describes the initial amplitude in r-space of the wave function. The number of elements is the number of internal degrees of freedom. Additional variables that are used in psi can be defined in the optional parameter variables. Predefined variables are the spatial coordinates x,y,z and time t.
  • number_of_grid_points (tuple of ints) – Tuple that defines the number of grid points (nX,nY,nZ) of the wave function
  • spatial_ext (list of tuples) – The supplied values define the boundary positions of the grid and thus define the actual coordinate system.
  • t0 (float, optional) – Internal time of wave function. Default is 0.0.
  • m (float, optional) – Mass of particle described by the wavefunction. Default is 1.054571817e-34 (numerically equal to hbar).
  • variables (dict) – Dictionary of additionally used variables in the definition of the wave function in psi. Predefined variables are the spatial coordinates x,y,z and time t.
  • normalize_const (float, optional) – Normalizes the wave function such that the integral of |Psi|^2 over all internal and external degrees of freedom equals normalize_const
num_int_dim

The number of internal degrees of freedom

Type:int
num_ex_dim

The number of external degrees of freedom

Type:int
r

The arrays are the evenly spaced spatial coordinates as defined through definition of spatial_ext and number_of_grid_points

Type:list of 1d arrays
k

The arrays are the evenly spaced inverse spatial coordinates as defined through definition of spatial_ext and number_of_grid_points

Type:list of 1d arrays

Examples

Wavefunction with two internal states where the first state is gaussian distributed in 1d r-space and the second state is not occupied at all.

>>> from pytalises import Wavefunction
>>> psi = Wavefunction(["exp(-((x-x0)/a0)**2)", "0.0"],
    (16,), [(-2,2),], variables={'a0':1/2, 'x0':0})
>>> print(psi.num_int_dim)
2
>>> print(psi.num_ext_dim)
1
>>> print(psi.amp)
[[1.12535175e-07+0.j 0.00000000e+00+0.j]
[6.03594712e-06+0.j 0.00000000e+00+0.j]
[1.83289361e-04+0.j 0.00000000e+00+0.j]
[3.15111160e-03+0.j 0.00000000e+00+0.j]
[3.06707930e-02+0.j 0.00000000e+00+0.j]
[1.69013315e-01+0.j 0.00000000e+00+0.j]
[5.27292424e-01+0.j 0.00000000e+00+0.j]
[9.31358402e-01+0.j 0.00000000e+00+0.j]
[9.31358402e-01+0.j 0.00000000e+00+0.j]
[5.27292424e-01+0.j 0.00000000e+00+0.j]
[1.69013315e-01+0.j 0.00000000e+00+0.j]
[3.06707930e-02+0.j 0.00000000e+00+0.j]
[3.15111160e-03+0.j 0.00000000e+00+0.j]
[1.83289361e-04+0.j 0.00000000e+00+0.j]
[6.03594712e-06+0.j 0.00000000e+00+0.j]
[1.12535175e-07+0.j 0.00000000e+00+0.j]]
>>> print(psi.r)
[array([-2.        , -1.73333333, -1.46666667, -1.2       , -0.93333333,
   -0.66666667, -0.4       , -0.13333333,  0.13333333,  0.4       ,
    0.66666667,  0.93333333,  1.2       ,  1.46666667,  1.73333333,
    2.        ]), array([0.]), array([0.])]
amp

Ndarray of the wave function amplitudes.

construct_FFT(num_of_threads=1, FFTWflags=('FFTW_ESTIMATE', 'FFTW_DESTROY_INPUT'))

Construct pyfftw bindings.

exp_pos(axis=None)

Calculate the expected position on given axis.

Calculates the mean position of Psi on chosen axis. Axes 0,1,2 correspond to x,y,z. The other two axes are traced out. If no axis iv given returns array of mean position of all external degrees of freedom.

freely_propagate(num_time_steps, delta_t, num_of_threads=1, FFTWflags=('FFTW_ESTIMATE', 'FFTW_DESTROY_INPUT'))

Propagates the Wavefunction object in time with V=0.

Function that can propagate the wavefunction if no potential is present.

Parameters:
  • num_time_steps (int) – Number of times the wavefunction is propagated by time delta_t using the Split-Step Fourier method.
  • delta_t (float) – Time increment the wavefunction is propagated in one time step.
  • num_of_threads (int, optional) – Number of threads uses for calculation. Default is 1.
  • FFTWflags (tuple of strings) – Options for FFTW planning [1]. Default is (‘FFTW_ESTIMATE’, ‘FFTW_DESTROY_INPUT’,).

References

[1] http://www.fftw.org/fftw3_doc/Planner-Flags.html

k

(list of) array of the wave function position axes.

normalize_to(n_const)

Normalize the wave function.

Normalizes the wave function such that the integral of |Psi|^2 over all internal and external states equals n_const

propagate(potential, num_time_steps, delta_t, **kwargs)

Propagates the Wavefunction object in time.

Function that propagates the wavefunction using a Split-Step Fourier method [1].

Parameters:
  • potential (string list of strings) – This list contains the matrix elements of the potential term V in string format. If the potential has nondiagonal elements (see optional parameter diag) each elements represents one matrix element of the lower triangular part of V. For example a 3x3 potential with nondiagonal elements would be of form potential=[H00, H10, H20, H11, H21, H22]. If the potential term is supposed to have only diagonal elements (diag=True), the potential parameter for a 3x3 potential would look like potential=[H00,H11,H22].
  • num_time_steps (int) – Number of times the wavefunction is propagated by time delta_t using the Split-Step Fourier method.
  • delta_t (float) – Time increment the wavefunction is propagated in one time step.
  • variables (dict, optional) – Dictionary containing values for variables you might have used in potential
  • diag (bool , optional) – If true, no numerical diagonalization has to be invoked in order to calculate time-propagation as nondiagonal elements are omitted. This makes the computation much faster. Default is False.
  • num_of_threads (int, optional) – Number of threads uses for calculation. Default behaviour is to use all threads available.
  • FFTWflags (tuple of strings) – Options for FFTW planning [2]. Default is (‘FFTW_ESTIMATE’, ‘FFTW_DESTROY_INPUT’,).

References

[1] https://en.wikipedia.org/wiki/Split-step_method [2] http://www.fftw.org/fftw3_doc/Planner-Flags.html

r

(list of) array of the wave function position axes.

state_occupation(nth_state=None)

Return occupation number of nth internal state.

Evaluates the spatial integral over |Psi|^2 for the nth internal state. If none is given a vector of the occupation number of all internal states is returned.

var_pos(axis=None)

Calculate the variance on given axis.

Calculates the variance in position of Psi on chosen axis. Axes 0,1,2 correspond to x,y,z. The other two axes are traced out. If no axis iv given returns array of variance position of all external degrees of freedom.

pytalises.wavefunction.assert_type_or_list_of_type(argument, wished_type)

pytalises.propagator module

Module containing functions that help propagating the Wavefunction class.

class pytalises.propagator.Propagator(psi, potential, variables={}, diag=False, num_of_threads=1, FFTWflags=('FFTW_ESTIMATE', 'FFTW_DESTROY_INPUT'))

Bases: object

Class for propagating instances of the Wavefunction class.

Parameters:
  • psi (Wavefunction) – The Wavefunction object the Propagator class acts on
  • potential (list of strings) – This list contains the matrix elements of the potential term V in string format. If the potential has nondiagonal elements (see optional parameter diag) each elements represents one matrix element of the lower triangular part of V. For example a 3x3 potential with nondiagonal elements would be of form potential=[H00, H10, H20, H11, H21, H22]. If the potential term is supposed to have only diagonal elements (diag=True), the potential argument for a 3x3 potential would look like potential=[H00,H11,H22].
  • variables (dict, optional) – Dictionary containing values for variables you might have used in potential
  • diag (bool , optional) – If true, no numerical diagonalization has to be invoked in order to calculate time-propagation. Default is False.
  • num_of_threads (int, optional) – Number of threads uses for calculation. Default is 1.
  • FFTWflags (tuple of strings) – Options for FFTW planning [1]. Default is (‘FFTW_ESTIMATE’, ‘FFTW_DESTROY_INPUT’,).

References

[1] http://www.fftw.org/fftw3_doc/Planner-Flags.html

class Potential(potential_string, variables={}, diag=False)

Bases: object

Simple class for collecting information about the potential.

diag_potential_prop(delta_t)

Calculate exp(i*V/hbar*delta_t)*Psi by simple matrix multiplication.

This method is used if the potential matrix V is diagonal. This is much faster than nondiag_potential_prop and should be used if possible.

eval_V()

Evalutes V on the whole spatial grid.

The result is saved in Propagator.V_eval_array.

eval_diag_V()

Evalutes diagonal elements of V on the whole spatial grid.

The result is saved in Propagator.V_eval_array.

kinetic_prop(delta_t)

Perform time propagation in k-space.

Transforms the Wavefunction into k-space, calculates exp(i*hbar/(2m)*k**2*delta_t)*Psi(kx,ky,kz) and transforms it back into r-space.

nondiag_potential_prop(delta_t)

Calculate exp(i*V/hbar*delta_t)*Psi using numerical diagonalization.

This method has to be used if the potential matrix has nondiagonal elements.

potential_prop(delta_t)

Wrap function that calculates exp(i*V(x,y,z)/hbar*delta_t)*Psi(x,y,z).

This can be either nondiag_potential_prop or diag_potential_prop.

pytalises.propagator.freely_propagate(psi, num_time_steps, delta_t, num_of_threads=1, FFTWflags=('FFTW_ESTIMATE', 'FFTW_DESTROY_INPUT'))

Propagates a Wavefunction object in time with V=0.

Function that can propagate the wavefunction if no potential is present.

Parameters:
  • psi (Wavefunction) – The Wavefunction object the Propagator class acts on
  • num_time_steps (int) – Number of times the wavefunction is propagated by time delta_t using the Split-Step Fourier method.
  • delta_t (float) – Time increment the wavefunction is propagated in one time step.
  • num_of_threads (int, optional) – Number of threads uses for calculation. Default is 1.
  • FFTWflags (tuple of strings) – Options for FFTW planning [1]. Default is (‘FFTW_ESTIMATE’, ‘FFTW_DESTROY_INPUT’,).

References

[1] http://www.fftw.org/fftw3_doc/Planner-Flags.html

pytalises.propagator.get_eig

Calculate eigenvectors and eigenvalues of matrices in array.

JIT-compiled function that calculates the eigenvectors and eigenvalues of input array M in parallel using numba. The resulting eigenvectors are stored in the input matrix and the eigenvalues in the array eigvals.

Parameters:
  • M (3d array of (NxN) arrays) –
  • eigvals (3d array of 1d arrays with N elements) –
pytalises.propagator.propagate(psi, potential, num_time_steps, delta_t, **kwargs)

Propagates a Wavefunction object in time.

Function that propagates the wavefunction using a Split-Step Fourier method [1].

Parameters:
  • psi (Wavefunction) – The Wavefunction object the Propagator class acts on
  • potential (string list of strings) – This list contains the matrix elements of the potential term V in string format. If the potential has nondiagonal elements (see optional parameter diag) each elements represents one matrix element of the lower triangular part of V. For example a 3x3 potential with nondiagonal elements would be of form potential=[H00, H10, H20, H11, H21, H22]. If the potential term is supposed to have only diagonal elements (diag=True), the potential parameter for a 3x3 potential would look like potential=[H00,H11,H22].
  • num_time_steps (int) – Number of times the wavefunction is propagated by time delta_t using the Split-Step Fourier method.
  • delta_t (float) – Time increment the wavefunction is propagated in one time step.
  • variables (dict, optional) – Dictionary containing values for variables you might have used in potential
  • diag (bool , optional) – If true, no numerical diagonalization has to be invoked in order to calculate time-propagation as nondiagonal elements are omitted. This makes the computation much faster. Default is False.
  • num_of_threads (int, optional) – Number of threads uses for calculation. Default behaviour is to use all threads available.
  • FFTWflags (tuple of strings) – Options for FFTW planning [2]. Default is (‘FFTW_ESTIMATE’, ‘FFTW_DESTROY_INPUT’,).

References

[1] https://en.wikipedia.org/wiki/Split-step_method [2] http://www.fftw.org/fftw3_doc/Planner-Flags.html