hqs_distiller.chebyshev

Module providing ChebyshevIterator class.

This module defines the ChebysgevIterator class, which facilitates the efficient computation of Chebyshev moments for the HQS Distiller package.

@private

Copyright © 2022-2025 HQS Quantum Simulations GmbH. All Rights Reserved.

class ChebyshevIterator:

Simple class implementing efficient Chebyshev recursion.

Initialization requires a scipy sparse matrix. During initialization this matrix is converted (copied!) to CSR format. Its shape is denoted by (L, L), 'L' being the linear size of the matrix.

NOTE: The spectrum of the matrix is supposed to be in the interval (-1, 1), otherwise the Chebyshev recursion is not stable. The user is responsible for rescaling the spectrum!!!

The __call__ method of the object expects two dense numpy arrays which correspond to the (n-1)'th and (n-2)'th Chebyshev polynomial of the matrix applied to the input state(s) (in the following called "Chebyshev moment(s)"). It returns the n'th Chebyshev moment(s).

NOTE: The __call__ method overwrites the (n-2)'th Chebyshev moment(s) by the n'th Chebyshev moment(s).

The apply_M method expects a dense numpy array, representing state(s) to which the linear operator M is applied. It returns the state(s) after application of the linear operator M. For example, this method can be used to generate the 1'st Chebyshev moment(s) [0'th -> trivially the input state(s)].

NOTE: Both methods can be applied to a set of N states, i.e., input numpy arrays of shape (L, N).

ChebyshevIterator(matrix: scipy.sparse._matrix.spmatrix)

Initialize Chebyshev iterator.

NOTE: Linear operator M must have spectrum in (-1, 1). Not checked in initialization!

NOTE: Only hermitian matrices are supported. Not checked in initialization!

The linear operator is stored in scipy's CSR format.

Arguments:
  • matrix (SPM.spmatrix): Linear operator for Chebyshev iteration.
Raises:
  • ValueError: If matrix is not a square matrix.
def apply_op(self, b: numpy.ndarray) -> numpy.ndarray:

Apply matrix M to given state(s).

NOTE: Linear operator A has shape (L, L).

Helper method to generate 1'st Chebyshev moment(s).

Arguments:
  • b (np.ndarray): State(s) to which M is applied. Shape must be (L, N).
Returns:

np.ndarray: State(s) after application of matrix M.

class ChebyshevIteratorOperator:

Simple class implementing efficient Chebyshev recursion.

Initialization requires a linear operator implemented as a Callable[[xp.ndarray, xp.ndarray], xp.ndarray] function. The signature implies that the first argument is the state to which the operator is applied. The second argument holds the (pre-allocated) output state. The return value is the output state.

NOTE: We use the "xp Pattern" to write the code to be compatible with both CPU (numpy) and GPU (cupy) environments.

NOTE: The spectrum of the matrix is supposed to be in the interval (-1, 1), otherwise the Chebyshev recursion is not stable. The user is responsible for rescaling the spectrum!!!

The __call__ method of the object expects two dense numpy arrays which correspond to the (n-1)'th and (n-2)'th Chebyshev polynomial of the matrix applied to the input state(s) (in the following called "Chebyshev moment(s)"). It returns the n'th Chebyshev moment(s).

NOTE: The __call__ method overwrites the (n-2)'th Chebyshev moment(s) by the n'th Chebyshev moment(s).

The apply_M method expects a dense numpy array, representing state(s) to which the linear operator M is applied. It returns the state(s) after application of the linear operator M. For example, this method can be used to generate the 1'st Chebyshev moment(s) [0'th -> trivially the input state(s)].

NOTE: Both methods can be applied to a set of N states, i.e., input numpy arrays of shape (L, N).

ChebyshevIteratorOperator( lin_op: Callable[[numpy.ndarray, numpy.ndarray], numpy.ndarray], state_shape: tuple[int, int], batch_size: int = 1)

Initialize Chebyshev iterator.

NOTE: Linear operator must have spectrum in (-1, 1). Not checked in initialization!

NOTE: Only hermitian matrices are supported. Not checked in initialization!

Arguments:
  • lin_op (Callable[[xp.ndarray, xp.ndarray], xp.ndarray]): Linear operator for Chebyshev iteration.
  • state_shape (tuple[int, int,]): Shape of state compatible with linear operator.
  • batch_size (int): Number of states in batch. Defaults to 1.
def apply_op(self, b: numpy.ndarray) -> numpy.ndarray:

Apply linear operator to given state.

NOTE: Linear operator A has shape (L, L).

Helper method to generate 1'st Chebyshev moment.

Arguments:
  • b (np.ndarray): State to which the linear operator is applied. Shape must be state_shape.
Returns:

np.ndarray: State after application of linear operator.

def get_projector_coefficients( n: int, window: tuple[float, float], character: str = 'ret') -> numpy.ndarray:

Returns Chebyshev expansion coefficients for resolvent.

FIXME: Function has considerable overlap with _get_DOS.

Arguments:
  • n (int): Number of Chebyshev moments.
  • window (tuple[float, float]): Energy window in 'Chebyshev interval'.
  • character (str): Specify convention for real. Currently supported: 'ret' -> sign(Im[z] = 0) = 1, 'adv' -> sign(Im[z] = 0) = -1. Defaults to 'ret'.
Returns:

np.ndarray: Expansion coefficients for resolvent.

Raises:
  • ValueError: If 'type' argument is not supported.
def get_density_of_states( chebyshev_moments: numpy.ndarray, window: tuple[float, float], chebyshev_window: tuple[float, float], n_grid_points: int = 100) -> tuple[float, numpy.ndarray, numpy.ndarray]:

Returns the density of states (DOS) in given energy window.

NOTE: The input requires precomputed Chebyshev moments of the Hamiltonian for which the DOS is computed.

NOTE: At present, both, the energy window and the energy window scaled into the Chebyshev window are required input.

NOTE: If 'bare' Chebyshev moments are used the DOS is normalized to 1. You can use 'pre-scaled' Chebyshev moments to enrforce a different normalization condition.

Arguments:
  • chebyshev_moments (np.ndarray): Precomputed Chebyshev moments.
  • window (tuple[float, float]): Energy window for which to compute the DOS.
  • chebyshev_window (tuple[float, float]): Energy window scaled into Chebyshev interval.
  • n_grid_points (int): Number of discretization points for energy window. Defaults to 100.
Returns:

tuple[float, np.ndarray, np.ndarray]: Tuple of float: energy spacing between grid points, np.ndarray: discretized energy window (equidistant grid), np.ndarray: DOS on grid points.

def initialize_states( cmg: ChebyshevIterator, n_chebyshev_moments: int, chebyshev_window: tuple[float, float], n_states: int = 8, verbose: bool = True) -> tuple[numpy.ndarray, numpy.ndarray, numpy.ndarray]:

Initialize states for FEAST algorithm.

Arguments:
  • cmg (ChebyshevIterator): Generator for Chebyshev recursion.
  • n_chebyshev_moments (int): Number of Chebyshev moments to be used.
  • chebyshev_window (tuple[float, float]): Energy window in 'Chebyshev interval' for the projection.
  • n_states (int): Number of states to be initialized. Defaults to 8.
  • verbose (bool): Flag indicating whether to print debug information. Default is True.
Returns:

tuple[np.ndarray, np.ndarray, np.ndarray]: Tuple of array of projected states, Chebyshev moments, trace of projector

def initialize_state( cmg: ChebyshevIteratorOperator, n_chebyshev_moments: int, chebyshev_window: tuple[float, float], dtype: numpy.dtype, verbose: bool = True) -> tuple[numpy.ndarray, numpy.ndarray]:

Initialize state for HQS Distiller.

Arguments:
  • cmg (ChebyshevIteratorOperator): Generator for Chebyshev recursion.
  • n_chebyshev_moments (int): Number of Chebyshev moments to be used.
  • chebyshev_window (tuple[float, float]): Energy window in 'Chebyshev interval' for the projection.
  • dtype (xp.dtype): Data type for the random state.
  • verbose (bool): Flag indicating whether to print debug information. Default is True.
Returns:

tuple[xp.ndarray, xp.ndarray]: Tuple of projected states, trace estimate of projector

def initialize_moments( cmg: ChebyshevIteratorOperator, n_chebyshev_moments: int, current_moments: numpy.ndarray, initial_state_c: numpy.ndarray, v0: numpy.ndarray, v1: numpy.ndarray, dtype: numpy.dtype, verbose: bool = True) -> tuple[numpy.ndarray, numpy.ndarray, numpy.ndarray]:

Initialize Chebyshev moments.

Arguments:
  • cmg (ChebyshevIteratorOperator): Generator for Chebyshev recursion.
  • n_chebyshev_moments (int): Number of Chebyshev moments to be used.
  • current_moments (xp.ndarray): Moments calculated in previous Chebyshev iteration.
  • initial_state_c (xp.ndarray): (Conjugate of) initial state for Chebyshev iteration.
  • v0 (xp.ndarray): second to last state of previous Chebyshev iteration.
  • v1 (xp.ndarray): last state in previous Chebyshev iteration.
  • dtype (xp.dtype): Data type supported by ChebyshevIteratorOperator.
  • verbose (bool): Flag indicating whether to print debug information. Default is True.
Returns:

tuple[xp.ndarray, xp.ndarray, xp.ndarray]: Tuple of Chebyshev moments, next to last state of chebyshev iteration, last state of chebyshev iteration.

def initialize_state_with_moments( cmg: ChebyshevIteratorOperator, n_chebyshev_moments: int, chebyshev_window: tuple[float, float], dtype: numpy.dtype, verbose: bool = True) -> tuple[numpy.ndarray, numpy.ndarray, numpy.ndarray]:

Initialize state for FEAST algorithm.

Arguments:
  • cmg (ChebyshevIteratorOperator): Generator for Chebyshev recursion.
  • n_chebyshev_moments (int): Number of Chebyshev moments to be used.
  • chebyshev_window (tuple[float, float]): Energy window in 'Chebyshev interval' for the projection.
  • dtype (xp.dtype): Data type for the random state.
  • verbose (bool): Flag indicating whether to print debug information. Default is True.
Returns:

tuple[xp.ndarray, xp.ndarray, xp.ndarray]: Tuple of projected states, Chebyshev moments, trace estimate of projector