dynamite.computations

dynamite.computations.evolve(H, state, t, result=None, tol=None, ncv=None, algo=None, max_its=None)[source]

Evolve a quantum state according to the Schrodinger equation under the Hamiltonian H. The units are natural, that is, the evolution is simply

\[\Psi_t = e^{-iHt} \Psi_0\]
Parameters:
  • H (Operator) – The Hamiltonian

  • state (dynamite.states.State) – A dynamite State object containing the state to be evolved.

  • t (float) – The time for which to evolve. Can be negative to evolve backwards in time.

  • result (dynamite.states.State, optional) – Where to store the result state. If not given, a new vector is created in which to store the result. If evolving repeatedly many times, it is a good idea to pass a result vector to avoid repeatedly allocating a lot of memory. Will be overwritten.

  • tol (float, optional) – The tolerance for the evolution. Error estimation is difficult for Krylov exponentiation; this merely requests that the error be somewhat close to tol. There is no guarantee that it will actually be smaller.

  • ncv (int, optional) – The Krylov subspace size to use. Increasing subspace size can increase performance by reducing the number of iterations necessary, but also linearly increases memory usage and the number of matrix multiplies performed. Optimizing this parameter can significantly affect performance.

  • algo (string, optional) – Allowed options: ‘krylov’ or ‘expokit’. Which SLEPc algorithm to use to compute the matrix exponential. Default is ‘expokit’.

  • max_its (int, optional) – Maximum number of iterations for the solver.

Returns:

The result state

Return type:

dynamite.states.State

dynamite.computations.eigsolve(H, getvecs=False, nev=1, which='smallest', target=None, tol=None, subspace=None, max_its=None)[source]

Solve for a subset of the eigenpairs of the Hamiltonian.

By default, solves for the eigenvalue with the smallest (most negative) real part, e.g. the ground state. Which eigenvalues are sought and how many can be adjusted with the options.

Note

Krylov algorithms have difficulty with degenerate or very nearly degenerate eigenvalues. Degenerate eigenvalues may be missed, and near-degenerate eigenstates may be inaccurate.

Note

Do not try to use this function to solve for the whole spectrum! It’s very efficient at finding a few eigenvalues, but no faster than other routines for finding all of them. In the future an efficient solver for the whole spectrum may be included with dynamite.

Parameters:
  • getvecs (Bool) – Whether to return eigenvectors as well as eigenvalues.

  • nev (int) – The number of eigenvalues sought. The algorithm may return more eigenvalues than nev if more happen to converge.

  • which (str) –

    Which eigenvalues to seek. Options are:

    • "smallest", to find the eigenvalues with smallest real part (i.e. most negative)

    • "largest", to find the eigenvalues with largest real part (i.e. most positive)

    • "exterior", to find eigenvalues largest in absolute magnitude

    • "target", to find eigenvalues closest to the given target

    If target is set, which can be omitted and will automatically be set to "target".

  • target (float) – Using the shift-invert method, the eigensolver can seek the eigenvalues with real part closest to some target value. This requires a linear solve and so will be slower than solving for exterior eigenvalues. PETSc must be configured with a parallel linear solver (e.g. --download-mumps option in configure) to use this option in parallel.

  • tol (float) – The tolerance for the computation.

  • subspace (dynamite.subspaces.Subspace, optional) – The subspace on which to solve for eigenvalues. If not given, defaults to the most recent subspace set with Operator.add_subspace, or config.subspace if no subspaces have been added.

Returns:

Either a 1D numpy array of eigenvalues, or a pair containing that array and a list of the corresponding eigenvectors.

Return type:

numpy.array or tuple(numpy.array, list(dynamite.states.State))

dynamite.computations.reduced_density_matrix(state, keep)[source]

Compute the reduced density matrix of a state vector by tracing out some set of spins. The spins to be kept (not traced out) are specified in the keep array.

The density matrix is returned on process 0, the function returns a 1x1 matrix containing the value -1 on all other processes.

Parameters:
  • state (dynamite.states.State) – A dynamite State object.

  • keep (array-like) – A list of spin indices to keep. Must be sorted. Note that the returned matrix will have dimension \(2^{\mathrm{len(keep)}}\), so too long a list will generate a huge matrix.

Returns:

The density matrix

Return type:

numpy.ndarray[np.complex128]

dynamite.computations.entanglement_entropy(state, keep)[source]

Compute the entanglement of a state across some cut on the spin chain. To be precise, this is the bipartite entropy of entanglement.

Currently, this quantity is computed entirely on process 0. As a result, the function returns -1 on all other processes.

Parameters:
Returns:

The entanglement entropy

Return type:

float

dynamite.computations.dm_entanglement_entropy(dm)[source]

Compute the Von Neumann entropy of a density matrix.

Parameters:

dm (np.array) – A density matrix

Returns:

The Von Neumann entropy

Return type:

float

dynamite.computations.renyi_entropy(state, keep, alpha, method='eigsolve')[source]

Compute the Renyi entropy of the density matrix that results from tracing out some spins. The Renyi entropy is

\[H_{\alpha} = \frac{1}{1 - \alpha} \log \mathrm{Tr} \left[ \rho ^ \alpha \right]\]

Arbitrary non-negative values of alpha are allowed; in the special cases of \(\alpha \in \{ 0, 1 \}\) the function is computed in the limit.

Currently, this quantity is computed entirely on process 0. As a result, the function returns -1 on all other processes.

Parameters:
  • state (dynamite.states.State) – A dynamite State object.

  • keep (array-like) – A list of spin indices to keep. See reduced_density_matrix() for details.

  • alpha (float, int, or str) – The value of \(\alpha\) from the definition of Renyi entropy.

  • method (str, optional) – Whether to compute the Renyi entropy by solving for eigenvalues, or computing a matrix power and doing a trace. One or the other may be faster depending on the specific problem. Options: eigsolve or matrix_power.

Returns:

The Renyi entropy

Return type:

float

dynamite.computations.dm_renyi_entropy(dm, alpha, method='eigsolve')[source]

Compute the Renyi entropy of a density matrix. See renyi_entropy() for details.

Parameters:
  • dm (np.array) – A density matrix

  • alpha (int, float, or str) – The value of alpha in the definition of Renyi entropy.

  • method (str) – Whether to compute the Renyi entropy by solving for eigenvalues, or computing a matrix power and doing a trace. One or the other may be faster depending on the specific problem. Options: eigsolve or matrix_power.

Returns:

The Renyi entropy

Return type:

float

dynamite.computations.get_tstep(ncv, nrm, tol=1e-07)[source]

Compute the length of a sub-step in a Expokit matrix exponential solve.

dynamite.computations.estimate_compute_time(t, ncv, nrm, tol=1e-07)[source]

Estimate compute time in units of matrix multiplies, for an expokit exponential solve.

exception dynamite.computations.ConvergenceError[source]

Bases: Exception

exception dynamite.computations.MaxIterationsError[source]

Bases: ConvergenceError