QSP, QSVT, and Block Encoding

Block encoding embeds a matrix in a unitary. QSP and QSVT implement polynomial or scalar-function transformations on that representation. This page distinguishes low-level circuit construction from higher-level result objects.

Block encoding


import numpy as np
from unitarylab.library import block_encode
 
A = np.array([[0.5, 0.1], [0.2, 0.3]], dtype=complex)
encoded = block_encode(A, method="nagy")
 
print(encoded.circuit)
print(encoded.alpha)
print(encoded.total_qubits, encoded.target_qubits)


block_encode(matrix, method="nagy", eps=1e-3, verbose=False)

Method	Characteristics	Parameters
`nagy`	Dilation-based block encoding, suitable for general small matrices	`eps` is not used by this method
`fable`	Fast approximate block encoding that can compress small rotations	`eps` is the compression threshold; `None` disables compression

Input may be a two-dimensional array, list, or tuple; a one-dimensional input is treated as a column vector. BlockEncodingResult contains circuit, alpha, total_qubits, target_qubits, method, matrix, and eps, and provides get_encoded_matrix(), get_unitary_matrix(), and get_max_error(). Combine multiple results with unitarylab.library.block_encoding.chain_matmul().

Nagy Numerical Precision

Nagy first normalizes and pads the input matrix, then constructs a unitary dilation. The current implementation computes both defect matrices from the same singular value decomposition, $A = V Σ W^{†}$ :

D = W max (I - Σ^{2}, 0) W^{†}, D_{1} = V max (I - Σ^{2}, 0) V^{†}, U = (A D D_{1} - A^{†}) .

Here $A$ is the normalized, padded matrix. max clips diagonal entries corresponding to the singular values elementwise to avoid taking square roots of negative roundoff values. Sharing the singular values and square roots fixes a numerical consistency issue caused by constructing the two defect matrices separately, particularly relevant when singular values approach 1. Internal defect-matrix computation uses complex double precision; this is not a fixed error guarantee for arbitrary inputs, precisions, or backends.

The block_encode() call, normalization factor alpha, and public BlockEncodingResult fields are unchanged. eps still does not control Nagy error. For small matrices, check both unitarity of the full dilation and recovery of the input matrix:


U = encoded.get_unitary_matrix()
unitarity_error = np.max(np.abs(U.conj().T @ U - np.eye(U.shape[0])))
restored = encoded.get_encoded_matrix()
print(unitarity_error)
assert np.allclose(restored, A, atol=1e-10)

get_encoded_matrix() returns an approximation of the input matrix after multiplying back by alpha; the top-left block of the full unitary corresponds to the normalized matrix. unitarity_error and input-matrix recovery error are separate checks.

QSVT


import numpy as np
from unitarylab.library import QSVT
 
H = np.array([[0.8, 0.0], [0.0, 0.4]], dtype=complex)
result = QSVT(
    H, function=lambda x: x * x, target_error=1e-5,
    block_encoding_method="nagy", backend="torch", device="cpu",
)
 
print(result.degree)
print(result.circuit)
print(result.total_error)

The high-level signature is:


QSVT(H, function, target_error=1e-6,
     block_encoding_method="nagy",
     backend="torch", device="cpu", dtype=np.complex128)

H must be Hermitian. function must accept NumPy scalars or arrays, and target_error must be positive. The returned QSVTResult provides H, function, degree, circuit, total_error, and evolution_result; the final three are computed on demand. For control over fitting stages, use unitarylab.library._qsvt.algorithm.QSVTHermitian.

Low-level QSP


QSP(U, n, m, coef, parity,
    eps=1e-12, maxiter=200, is_coef_cheby=False)

U: signal unitary circuit to transform.
n: number of target-register qubits; m: polynomial-order-related parameter.
coef: polynomial coefficients; parity: target polynomial parity.
If is_coef_cheby=False, coef contains ordinary polynomial coefficients; if True, it contains Chebyshev coefficients.
eps and maxiter control phase-solving convergence.

The Hamiltonian-specific wrapper returns a five-element tuple:


circuit, factor, n_ancilla, n_qubits, degree = QSP_hamiltonian_simulation(
    U_H, n, alpha, m, t, epsilon, beta, flag
)

U_H is a block-encoding circuit for $H / α$ ; flag selects the construction branch. The input scale must agree with the block-encoding normalization factor alpha, or the transformation will not represent the target Hamiltonian.

Constraints and cost

Block-encoding matrices must be finite and two-dimensional; an unknown method raises ValueError.
Constructing the full unitary matrix and computing errors can substantially increase memory use.
QSP/QSVT degree grows with approximation precision and function shape; validate with low degree and small matrices first.
When composing block encodings, adjacent result matrix dimensions must be compatible.