quimb.linalg.approx_spectral

Use stochastic Lanczos quadrature to approximate spectral function sums of any operator which has an efficient representation of action on a vector.

Attributes

njit

Numba no-python jit, but obeying cache setting.

ptr

Alias for partial_trace().

reqs

Functions

divide_update_(X, c, out)

Accelerated computation of X / c into out.

dot(a, b)

Matrix multiplication, dispatched to dense or sparse functions.

prod(iterable)

subtract_update_(X, c, Y)

Accelerated inplace computation of X -= c * Y. This is mainly

vdot(a, b)

Accelerated 'Hermitian' inner product of two arrays. In other words,

rand_phase(shape[, scale, dtype])

Generate random complex numbers distributed on the unit sphere.

rand_rademacher(shape[, scale, loc, dtype])

randn([shape, dtype, scale, loc, num_threads, seed, dist])

Fast multithreaded generation of random normally distributed data.

seed_rand(seed)

See the random number generators, by instantiating a new set of bit

get_mpi_pool([num_workers, num_threads])

Get the MPI executor pool, with specified number of processes and

default_to_neutral_style(fn)

Wrap a function or method to use the neutral style by default.

find_library(x)

Check if library is installed.

format_number_with_error(x, err)

Given x with error err, format a string showing the relevant

int2tup(x)

raise_cant_find_library_function(x[, extra_msg])

Return function to flag up a missing necessary library.

lazy_ptr_linop(psi_ab, dims, sysa, **linop_opts)

A linear operator representing action of partially tracing a bipartite

lazy_ptr_ppt_linop(psi_abc, dims, sysa, sysb, **linop_opts)

A linear operator representing action of partially tracing a tripartite

inner(a, b)

Inner product between two vectors

norm_fro(a)

'Frobenius' norm of a vector.

norm_fro_approx(A, **kwargs)

Calculate the approximate frobenius norm of any hermitian linear

random_rect(shape[, dist, orthog, norm, seed, dtype])

Generate a random array optionally orthogonal.

construct_lanczos_tridiag(A, K[, v0, bsz, k_min, ...])

Construct the tridiagonal lanczos matrix using only matvec operators.

lanczos_tridiag_eig(alpha, beta[, check_finite])

Find the eigen-values and -vectors of the Lanczos triadiagonal matrix.

calc_trace_fn_tridiag(tl, tv, f[, pos])

Spectral ritz function sum, weighted by ritz vectors.

ext_per_trim(x[, p, s])

Extended percentile trimmed-mean. Makes the mean robust to asymmetric

nbsum(xs)

std(xs)

Simple standard deviation - don't invoke numpy for small lists.

calc_est_fit(estimates, conv_n, tau)

Make estimate by fitting exponential convergence to estimates.

calc_est_window(estimates, conv_n)

Make estimate from mean of last m samples, following:

single_random_estimate(A, K, bsz, beta_tol, v0, f, ...)

calc_stats(samples, mean_p, mean_s, tol, tol_scale)

Get an estimate from samples.

get_single_precision_dtype(dtype)

get_equivalent_real_dtype(dtype)

plot_approx_spectral_info(info)

approx_spectral_function(A, f[, tol, bsz, R, R_min, ...])

Approximate a spectral function, that is, the quantity Tr(f(A)).

tr_abs_approx(*args, **kwargs)

tr_exp_approx(*args, **kwargs)

tr_sqrt_approx(*args, **kwargs)

xlogx(x)

tr_xlogx_approx(*args, **kwargs)

entropy_subsys_approx(psi_ab, dims, sysa[, backend])

Approximate the (Von Neumann) entropy of a pure state's subsystem.

tr_sqrt_subsys_approx(psi_ab, dims, sysa[, backend])

Approximate the trace sqrt of a pure state's subsystem.

norm_ppt_subsys_approx(psi_abc, dims, sysa, sysb[, ...])

Estimate the norm of the partial transpose of a pure state's subsystem.

logneg_subsys_approx(psi_abc, dims, sysa, sysb, **kwargs)

Estimate the logarithmic negativity of a pure state's subsystem.

negativity_subsys_approx(psi_abc, dims, sysa, sysb, ...)

Estimate the negativity of a pure state's subsystem.

gen_bipartite_spectral_fn(exact_fn, approx_fn, ...)

Generate a function that computes a spectral quantity of the subsystem

Module Contents

quimb.linalg.approx_spectral.divide_update_(X, c, out)[source]

Accelerated computation of X / c into out.

quimb.linalg.approx_spectral.dot(a, b)[source]

Matrix multiplication, dispatched to dense or sparse functions.

Parameters:
  • a (dense or sparse operator) – First array.

  • b (dense or sparse operator) – Second array.

Returns:

Dot product of a and b.

Return type:

dense or sparse operator

quimb.linalg.approx_spectral.njit

Numba no-python jit, but obeying cache setting.

quimb.linalg.approx_spectral.prod(iterable)
quimb.linalg.approx_spectral.ptr[source]

Alias for partial_trace().

quimb.linalg.approx_spectral.subtract_update_(X, c, Y)[source]

Accelerated inplace computation of X -= c * Y. This is mainly for Lanczos iteration.

quimb.linalg.approx_spectral.vdot(a, b)[source]

Accelerated ‘Hermitian’ inner product of two arrays. In other words, b here will be conjugated by the function.

quimb.linalg.approx_spectral.rand_phase(shape, scale=1, dtype=complex)[source]

Generate random complex numbers distributed on the unit sphere.

quimb.linalg.approx_spectral.rand_rademacher(shape, scale=1, loc=0.0, dtype=float)[source]
quimb.linalg.approx_spectral.randn(shape=(), dtype=float, scale=1.0, loc=0.0, num_threads=None, seed=None, dist='normal')[source]

Fast multithreaded generation of random normally distributed data.

Parameters:
  • shape (tuple[int]) – The shape of the output random array.

  • dtype ({'complex128', 'float64', 'complex64' 'float32'}, optional) – The data-type of the output array.

  • scale (float, optional) – A multiplicative scale for the random numbers.

  • loc (float, optional) – An additive location for the random numbers.

  • num_threads (int, optional) – How many threads to use. If None, decide automatically.

  • dist ({'normal', 'uniform', 'rademacher', 'exp'}, optional) – Type of random number to generate.

quimb.linalg.approx_spectral.seed_rand(seed)[source]

See the random number generators, by instantiating a new set of bit generators with a ‘seed sequence’.

quimb.linalg.approx_spectral.get_mpi_pool(num_workers=None, num_threads=1)

Get the MPI executor pool, with specified number of processes and threads per process.

quimb.linalg.approx_spectral.default_to_neutral_style(fn)[source]

Wrap a function or method to use the neutral style by default.

quimb.linalg.approx_spectral.find_library(x)[source]

Check if library is installed.

Parameters:

x (str) – Name of library

Returns:

If library is available.

Return type:

bool

quimb.linalg.approx_spectral.format_number_with_error(x, err)[source]

Given x with error err, format a string showing the relevant digits of x with two significant digits of the error bracketed, and overall exponent if necessary.

Parameters:
  • x (float) – The value to print.

  • err (float) – The error on x.

Return type:

str

Examples

>>> print_number_with_uncertainty(0.1542412, 0.0626653)
'0.154(63)'
>>> print_number_with_uncertainty(-128124123097, 6424)
'-1.281241231(64)e+11'
quimb.linalg.approx_spectral.int2tup(x)[source]
quimb.linalg.approx_spectral.raise_cant_find_library_function(x, extra_msg=None)[source]

Return function to flag up a missing necessary library.

This is simplify the task of flagging optional dependencies only at the point at which they are needed, and not earlier.

Parameters:
  • x (str) – Name of library

  • extra_msg (str, optional) – Make the function print this message as well, for additional information.

Returns:

A mock function that when called, raises an import error specifying the required library.

Return type:

callable

quimb.linalg.approx_spectral.reqs = '[cotengra,autoray]'
quimb.linalg.approx_spectral.lazy_ptr_linop(psi_ab, dims, sysa, **linop_opts)[source]

A linear operator representing action of partially tracing a bipartite state, then multiplying another ‘unipartite’ state:

  ( | )
+-------+
| psi_a |   ______
+_______+  /      \
   a|      |b     |
+-------------+   |
|  psi_ab.H   |   |
+_____________+   |
                  |
+-------------+   |
|   psi_ab    |   |
+_____________+   |
   a|      |b     |
    |      \______/
Parameters:
  • psi_ab (ket) – State to partially trace and dot with another ket, with size prod(dims).

  • dims (sequence of int, optional) – The sub dimensions of psi_ab.

  • sysa (int or sequence of int, optional) – Index(es) of the ‘a’ subsystem(s) to keep.

quimb.linalg.approx_spectral.lazy_ptr_ppt_linop(psi_abc, dims, sysa, sysb, **linop_opts)[source]

A linear operator representing action of partially tracing a tripartite state, partially transposing the remaining bipartite state, then multiplying another bipartite state:

     ( | )
+--------------+
|   psi_ab     |
+______________+  _____
 a|  ____   b|   /     \
  | /   a\   |   |c    |
  | | +-------------+  |
  | | |  psi_abc.H  |  |
  \ / +-------------+  |
   X                   |
  / \ +-------------+  |
  | | |   psi_abc   |  |
  | | +-------------+  |
  | \____/a  |b  |c    |
 a|          |   \_____/
Parameters:
  • psi_abc (ket) – State to partially trace, partially transpose, then dot with another ket, with size prod(dims). prod(dims[sysa] + dims[sysb]).

  • dims (sequence of int) – The sub dimensions of psi_abc.

  • sysa (int or sequence of int, optional) – Index(es) of the ‘a’ subsystem(s) to keep, with respect to all the dimensions, dims, (i.e. pre-partial trace).

  • sysa – Index(es) of the ‘b’ subsystem(s) to keep, with respect to all the dimensions, dims, (i.e. pre-partial trace).

quimb.linalg.approx_spectral.inner(a, b)[source]

Inner product between two vectors

quimb.linalg.approx_spectral.norm_fro(a)[source]

‘Frobenius’ norm of a vector.

quimb.linalg.approx_spectral.norm_fro_approx(A, **kwargs)[source]

Calculate the approximate frobenius norm of any hermitian linear operator:

\[\mathrm{Tr} \left[ A^{\dagger} A \right]\]
Parameters:
  • A (linear operator like) – Operator with a dot method, assumed to be hermitian, to estimate the frobenius norm of.

  • kwargs – Supplied to approx_spectral_function().

Return type:

float

quimb.linalg.approx_spectral.random_rect(shape, dist='rademacher', orthog=False, norm=True, seed=False, dtype=complex)[source]

Generate a random array optionally orthogonal.

Parameters:
  • shape (tuple of int) – The shape of array.

  • dist ({'guassian', 'rademacher'}) – Distribution of the random variables.

  • orthog (bool or operator.) – Orthogonalize the columns if more than one.

  • norm (bool) – Explicitly normalize the frobenius norm to 1.

quimb.linalg.approx_spectral.construct_lanczos_tridiag(A, K, v0=None, bsz=1, k_min=10, orthog=False, beta_tol=1e-06, seed=False, v0_opts=None)[source]

Construct the tridiagonal lanczos matrix using only matvec operators. This is a generator that iteratively yields the alpha and beta digaonals at each step.

Parameters:
  • A (dense array, sparse matrix or linear operator) – The operator to approximate, must implement .dot method to compute its action on a vector.

  • K (int, optional) – The maximum number of iterations and thus rank of the matrix to find.

  • v0 (vector, optional) – The starting vector to iterate with, default to random.

  • bsz (int, optional) – The block size (number of columns) of random vectors to iterate with.

  • k_min (int, optional) – The minimum size of the krylov subspace for form.

  • orthog (bool, optional) – If True, perform full re-orthogonalization for each new vector.

  • beta_tol (float, optional) – The ‘breakdown’ tolerance. If the next beta ceofficient in the lanczos matrix is less that this, implying that the full non-null space has been found, terminate early.

  • seed (bool, optional) – If True, seed the numpy random generator with a system random int.

Yields:
  • alpha (sequence of float of length k) – The diagonal entries of the lanczos matrix.

  • beta (sequence of float of length k) – The off-diagonal entries of the lanczos matrix, with the last entry the ‘look’ forward value.

  • scaling (float) – How to scale the overall weights.

quimb.linalg.approx_spectral.lanczos_tridiag_eig(alpha, beta, check_finite=True)[source]

Find the eigen-values and -vectors of the Lanczos triadiagonal matrix.

Parameters:
  • alpha (array of float) – The diagonal.

  • beta (array of float) – The k={-1, 1} off-diagonal. Only first len(alpha) - 1 entries used.

quimb.linalg.approx_spectral.calc_trace_fn_tridiag(tl, tv, f, pos=True)[source]

Spectral ritz function sum, weighted by ritz vectors.

quimb.linalg.approx_spectral.ext_per_trim(x, p=0.6, s=1.0)[source]

Extended percentile trimmed-mean. Makes the mean robust to asymmetric outliers, while using all data when it is nicely clustered. This can be visualized roughly as:

    |--------|=========|--------|
x     x   xx xx xxxxx xxx   xx     x      x           x

Where the inner range contains the central p proportion of the data, and the outer ranges entends this by a factor of s either side.

Parameters:
  • x (array) – Data to trim.

  • p (Proportion of data used to define the 'central' percentile.) – For example, p=0.5 gives the inter-quartile range.

  • s (Include data up to this factor times the central 'percentile' range) – away from the central percentile itself.

  • Returns

  • xt (array) – Trimmed data.

quimb.linalg.approx_spectral.nbsum(xs)[source]
quimb.linalg.approx_spectral.std(xs)[source]

Simple standard deviation - don’t invoke numpy for small lists.

quimb.linalg.approx_spectral.calc_est_fit(estimates, conv_n, tau)[source]

Make estimate by fitting exponential convergence to estimates.

quimb.linalg.approx_spectral.calc_est_window(estimates, conv_n)[source]

Make estimate from mean of last m samples, following:

  1. Take between conv_n and 12 estimates.

  2. Pair the estimates as they are alternate upper/lower bounds

  3. Compute the standard error on the paired estimates.

quimb.linalg.approx_spectral.single_random_estimate(A, K, bsz, beta_tol, v0, f, pos, tau, tol_scale, k_min=10, verbosity=0, *, seed=None, v0_opts=None, info=None, **lanczos_opts)[source]
quimb.linalg.approx_spectral.calc_stats(samples, mean_p, mean_s, tol, tol_scale)[source]

Get an estimate from samples.

quimb.linalg.approx_spectral.get_single_precision_dtype(dtype)[source]
quimb.linalg.approx_spectral.get_equivalent_real_dtype(dtype)[source]
quimb.linalg.approx_spectral.plot_approx_spectral_info(info)[source]
quimb.linalg.approx_spectral.approx_spectral_function(A, f, tol=0.01, *, bsz=1, R=1024, R_min=3, tol_scale=1, tau=0.0001, k_min=10, k_max=512, beta_tol=1e-06, mpi=False, mean_p=0.7, mean_s=1.0, pos=False, v0=None, verbosity=0, single_precision='AUTO', info=None, progbar=False, plot=False, **lanczos_opts)[source]

Approximate a spectral function, that is, the quantity Tr(f(A)).

Parameters:
  • A (dense array, sparse matrix or LinearOperator) – Operator to approximate spectral function for. Should implement A.dot(vec).

  • f (callable) – Scalar function with which to act on approximate eigenvalues.

  • tol (float, optional) – Relative convergence tolerance threshold for error on mean of repeats. This can pretty much be relied on as the overall accuracy. See also tol_scale and tau. Default: 1%.

  • bsz (int, optional) – Number of simultenous vector columns to use at once, 1 equating to the standard lanczos method. If bsz > 1 then A must implement matrix-matrix multiplication. This is a more performant way of essentially increasing R, at the cost of more memory. Default: 1.

  • R (int, optional) – The number of repeats with different initial random vectors to perform. Increasing this should increase accuracy as sqrt(R). Cost of algorithm thus scales linearly with R. If tol is non-zero, this is the maximum number of repeats.

  • R_min (int, optional) – The minimum number of repeats to perform. Default: 3.

  • tau (float, optional) – The relative tolerance required for a single lanczos run to converge. This needs to be small enough that each estimate with a single random vector produces an unbiased sample of the operators spectrum..

  • k_min (int, optional) – The minimum size of the krylov subspace to form for each sample.

  • k_max (int, optional) – The maximum size of the kyrlov space to form. Cost of algorithm scales linearly with K. If tau is non-zero, this is the maximum size matrix to form.

  • tol_scale (float, optional) – This sets the overall expected scale of each estimate, so that an absolute tolerance can be used for values near zero. Default: 1.

  • beta_tol (float, optional) – The ‘breakdown’ tolerance. If the next beta ceofficient in the lanczos matrix is less that this, implying that the full non-null space has been found, terminate early. Default: 1e-6.

  • mpi (bool, optional) – Whether to parallelize repeat runs over MPI processes.

  • mean_p (float, optional) – Factor for robustly finding mean and err of repeat estimates, see ext_per_trim().

  • mean_s (float, optional) – Factor for robustly finding mean and err of repeat estimates, see ext_per_trim().

  • v0 (vector, or callable) – Initial vector to iterate with, sets R=1 if given. If callable, the function to produce a random intial vector (sequence).

  • pos (bool, optional) – If True, make sure any approximate eigenvalues are positive by clipping below 0.

  • verbosity ({0, 1, 2}, optional) – How much information to print while computing.

  • single_precision ({'AUTO', False, True}, optional) – Try and convert the operator to single precision. This can lead to much faster operation, especially if a GPU is available. Additionally, double precision is not really needed given the stochastic nature of the algorithm.

  • lanczos_opts – Supplied to single_random_estimate() or construct_lanczos_tridiag().

Returns:

The approximate value Tr(f(a)).

Return type:

scalar

quimb.linalg.approx_spectral.tr_abs_approx(*args, **kwargs)[source]
quimb.linalg.approx_spectral.tr_exp_approx(*args, **kwargs)[source]
quimb.linalg.approx_spectral.tr_sqrt_approx(*args, **kwargs)[source]
quimb.linalg.approx_spectral.xlogx(x)[source]
quimb.linalg.approx_spectral.tr_xlogx_approx(*args, **kwargs)[source]
quimb.linalg.approx_spectral.entropy_subsys_approx(psi_ab, dims, sysa, backend=None, **kwargs)[source]

Approximate the (Von Neumann) entropy of a pure state’s subsystem.

Parameters:
  • psi_ab (ket) – Bipartite state to partially trace and find entopy of.

  • dims (sequence of int, optional) – The sub dimensions of psi_ab.

  • sysa (int or sequence of int, optional) – Index(es) of the ‘a’ subsystem(s) to keep.

  • kwargs – Supplied to approx_spectral_function().

quimb.linalg.approx_spectral.tr_sqrt_subsys_approx(psi_ab, dims, sysa, backend=None, **kwargs)[source]

Approximate the trace sqrt of a pure state’s subsystem.

Parameters:
  • psi_ab (ket) – Bipartite state to partially trace and find trace sqrt of.

  • dims (sequence of int, optional) – The sub dimensions of psi_ab.

  • sysa (int or sequence of int, optional) – Index(es) of the ‘a’ subsystem(s) to keep.

  • kwargs – Supplied to approx_spectral_function().

quimb.linalg.approx_spectral.norm_ppt_subsys_approx(psi_abc, dims, sysa, sysb, backend=None, **kwargs)[source]

Estimate the norm of the partial transpose of a pure state’s subsystem.

quimb.linalg.approx_spectral.logneg_subsys_approx(psi_abc, dims, sysa, sysb, **kwargs)[source]

Estimate the logarithmic negativity of a pure state’s subsystem.

Parameters:
  • psi_abc (ket) – Pure tripartite state, for which estimate the entanglement between ‘a’ and ‘b’.

  • dims (sequence of int) – The sub dimensions of psi_abc.

  • sysa (int or sequence of int, optional) – Index(es) of the ‘a’ subsystem(s) to keep, with respect to all the dimensions, dims, (i.e. pre-partial trace).

  • sysa – Index(es) of the ‘b’ subsystem(s) to keep, with respect to all the dimensions, dims, (i.e. pre-partial trace).

  • kwargs – Supplied to approx_spectral_function().

quimb.linalg.approx_spectral.negativity_subsys_approx(psi_abc, dims, sysa, sysb, **kwargs)[source]

Estimate the negativity of a pure state’s subsystem.

Parameters:
  • psi_abc (ket) – Pure tripartite state, for which estimate the entanglement between ‘a’ and ‘b’.

  • dims (sequence of int) – The sub dimensions of psi_abc.

  • sysa (int or sequence of int, optional) – Index(es) of the ‘a’ subsystem(s) to keep, with respect to all the dimensions, dims, (i.e. pre-partial trace).

  • sysa – Index(es) of the ‘b’ subsystem(s) to keep, with respect to all the dimensions, dims, (i.e. pre-partial trace).

  • kwargs – Supplied to approx_spectral_function().

quimb.linalg.approx_spectral.gen_bipartite_spectral_fn(exact_fn, approx_fn, pure_default)[source]

Generate a function that computes a spectral quantity of the subsystem of a pure state. Automatically computes for the smaller subsystem, or switches to the approximate method for large subsystems.

Parameters:
  • exact_fn (callable) – The function that computes the quantity on a density matrix, with signature: exact_fn(rho_a, rank=...).

  • approx_fn (callable) – The function that approximately computes the quantity using a lazy representation of the whole system. With signature approx_fn(psi_ab, dims, sysa, **approx_opts).

  • pure_default (float) – The default value when the whole state is the subsystem.

Returns:

bipartite_spectral_fn – The function, with signature: (psi_ab, dims, sysa, approx_thresh=2**13, **approx_opts)

Return type:

callable