Spin-½ Operators

Operators are a central tool in expressing and solving quantum many body physics problems in Aleph. They underpin Aleph's ability to mix symbolic computation and numerical methods. In this document we will review the operator framework from the perspective of spin-½ systems. Operators make writing expressions of quantum objects simple, while avoiding costly memory allocations. In this document we will explore the full suite of operators available to users.

Purely symbolic operators

The first set of operators we care about are purely "symbolic" operators, that is, operators which do not need anything other than a symbol to describe them.

Operator (example)	Support size	Description
ID()	N/A	Identity operator
X(0)	1	Pauli-X $\sigma^x$
Y(0)	1	Pauli-Y $\sigma^y$
Z(0)	1	Pauli-Z $\sigma^z$
XX(0,1)	2	Two-qubit XX interaction
YY(0,1)	2	Two-qubit YY interaction
ZZ(0,1)	2	Two-qubit ZZ interaction
XXPYY	2	$XX(i,j) + YY(i,j)$ optimization
SX(0)	1	$\frac{1}{2}\sigma^x$
SY(0)	1	$\frac{1}{2}\sigma^y$
SZ(0)	1	$\frac{1}{2}\sigma^z$
Splus(0)	1	Spin-up (S⁺)
Sminus(0)	1	Spin-down (S⁻)
Proj0(0)	1	projector $\|0 \rangle \langle 0\|$
Proj1(0)	1	projector $\|1\rangle\langle 1\|$
Hgate(0)	1	Hadamard gate
Sgate(0)	1	S gate
Tgate(0)	1	T gate
SWAP(0,1)	2	SWAP gate
CNOT(0,1)	2	CNOT gate
ZZNN([0,1,2,4])	N	$ZZ$ nearest neighbour sum $\sum_i Z_i Z_{i+1}$
ZZNNN([0,1,2,4])	N	$ZZ$ next nearest neighbour sum $\sum_i Z_i Z_{i+2}$

Parameterized operators

The next class are special operators which require the user to specify numbers to fully characterize the operator. A number of such operators are supported.

Operator (example)	Support size	param type	Description
Phase($\phi$)	0	(real)	$e^{-i\phi}$
RotX($\theta$,0)	1	(real, integer)	$e^{-i\frac{\theta}{2}\sigma^x}$
RotY($\theta$,0)	1	(real, integer)	$e^{-i\frac{\theta}{2}\sigma^y}$
RotZ($\theta$,0)	1	(real, integer)	$e^{-i\frac{\theta}{2}\sigma^z}$
RotEuler($\omega$, $\theta$, $\phi$,0)	1	(real, real, real, integer)	Euler rotation
DualUnitary($J$, [0,1])	2	(real,[integer,integer])	$\exp\left(i \left[ \frac{\pi}{4} XX + \frac{\pi}{4} YY + JZZ \right] \right)$
UnitaryXYZ$(J_x, J_y, J_z)$	2	([real,real,real],[integer,integer])	$\exp(-i[J_x X X + J_y Y Y + J_z Z Z])$
Floquet([$J_x,J_y,J_z$],[0,1])	2	([real,real,real],[integer,integer])	$\exp \left( -i [J_x X + J_y Y + J_z Z ] Z \right)$
pauli_string("XYZX", [0,1,2,3])	N	(string, [integer,..])	String of Pauli matrices.
operator_diagonal(vector,[1,2,3], as_complex)	N	(matrix, [integer, ..], optional type info)	Diagonal operator represented by a vector
operator_dense(matrix,[0,1])	N	(matrix, [integer,..])	Local operator represented as a $2^N \times 2^N$ matrix
operator_sparse(matrix,[0,1])	N	(matrix, [integer,..])	Sparse operator represented by a matrix
Flip(0)	1	(integer)	Flips spins at a sites
FlipAll()	0	()	Flips all spins
Reflect()	0	()	Reflects all spins about the center
Translate(1)	1	(integer)	Translates spins by a specified shift

Combining operators into sums and products

Individual operators are useful, but quantum Hamiltonians, circuits and more complicated expressions require combinations of operators through products and sums. We provide two factory functions to build these combinations: operator_prod() for products and operator_sum() for sums.

Building operator products with operator_prod()

Use operator_prod() to create products of operators. You can initialize it with an expression or start empty. Products have no size limit and can be indexed to access individual operators.

// Create from an expression
var prod = operator_prod(X(0)*X(1)); // -> X(0) * X(1) 
prod *= Y(5); // -> X(0) * X(1) * Y(5) 

// Or start empty (behaves like identity)
var prod2 = operator_prod()
for(m : [0..1000])
{
    prod2 *= ZZ(m, m+1)
}
print(prod2[0]) // -> Z(0,1)

warning

operator_prod() does not take any coefficients. Please use operator_sum() instead if you want to scale an operator.

Building operator sums with operator_sum()

Use operator_sum() to create sums of operators with either strictly real or complex coefficients. By default, operator_sum() uses complex coefficients unless you explicitly specify as_real:

// Default: complex coefficients
var default_sum = operator_sum();  // Same as operator_sum(as_complex)

// Explicitly specify real coefficients
var realopsum = operator_sum(as_real);

// Explicitly specify complex coefficients
var complexopsum = operator_sum(as_complex);

// Build up the sums
realopsum += 2.0 * X(0)*Z(1) + 3.5*X(1);
complexopsum += (1.0 + 2.0i)*Y(2) + 0.5i*Z(3);

You can also initialize directly from an expression, and you may specify the coefficient type:

// Initialize with real coefficients
var H_real = operator_sum(2.0*X(0) + 3.0*Y(1), as_real);

// Initialize with complex coefficients  
var H_complex = operator_sum((1.0+1.0i)*X(0) + 2.0i*Y(1), as_complex);

Choosing coefficient types

Use as_real when all coefficients are real for better performance and memory efficiency.

Use as_complex (or omit the type specifier for default behavior) when you need complex coefficients.

Do not mix real and complex coefficients

When using complex coefficients, all coefficients in the expression are treated as complex, even if the imaginary part is zero.

// ❌ INCORRECT: Mixing real and complex coefficients
var H_bad = operator_sum(2.0 * X(0) + 3.0i * Y(1));
var H_still_bad = operator_sum(2.0 * X(0) + 3.0i * Y(1), as_complex);

// ✅ CORRECT: All coefficients in complex form
var H_good = operator_sum((2.0 + 0i) * X(0) + (0.0 + 3.0i) * Y(1));

For purely real coefficients, use as_real to avoid this requirement and gain better performance.

Coefficient management

Operator sums do not automatically merge like terms. For example, X(0) + X(0) stores two separate terms with coefficient 1.0 each. Use the merged() function (see Operator Transformations) to combine like terms when needed.

Accessing individual operators

You can access individual operators within products and sums using the [] accessor. Indices read from left to right (index 0 is the leftmost operator):

var prod = Z(0) * RotZ(pi/3, 1) * X(2);
print(prod[0]);  // -> Z(0)
print(prod[1]);  // -> RotZ(pi/3, 1)
print(prod[2]);  // -> X(2)

Similarly, for operator_sum, you can access individual terms:

var sum = operator_sum(2.0 * X(0) + 3.5 * Y(1) + 1.0 * Z(2), as_real);
print(sum[0].second);  // -> X(0)
print(sum[1].second);  // -> Y(1)
print(sum[2].second);  // -> Z(2)

Extracting coefficients in an operator_sum

You can extract the coefficients from an operator_sum using get_coefficients():

// Real operator sum
var H_real = operator_sum(as_real);
H_real += 2.0 * X(0) + 3.5 * Y(1);
var real_coeffs = get_coefficients(H_real);  // Returns [2, 3.5]
var rc0 = H_real[0].first; // 2.0
var rc1 = H_real[1].first; // 3.5

// Complex operator sum
var H_complex = operator_sum(as_complex);
H_complex += (1.0 + 2.0i) * X(0) + (3.0 - 1.5i) * Y(1);
var complex_coeffs = get_coefficients(H_complex);  // Returns [1 + 2i, 3 - 1.5i]
var cc0 = H_complex[0].first; // (1.0 + 2.0i)
var cc1 = H_complex[1].first; // (3.0 - 1.5i)

This is useful for analyzing the structure of your Hamiltonian or inspecting individual terms before applying transformations.

Getting operator sites

Use support() to get the set of sites that an operator acts on:

var op1 = X(0);
var sup1 = support(op1);  // Returns [0]

var op2 = X(0) * Y(2) * Z(5);
var sup2 = support(op2);  // Returns [5, 2, 0]

var op3 = CNOT(3, 7);
var sup3 = op3.support();  // Returns [7, 3]

Current limitation

This version supports operator_prod and single operators. Support for operator_sum may not work correctly in the current version.

Use sites() to get the sites of an individual Spin-½ operator in definition order:

var op1 = X(0);
var sites1 = sites(op1);  // Returns [0]

var op2 = XX(0, 1);
var sites2 = sites(op2);  // Returns [0, 1]

var op3 = CNOT(2, 0);
var sites3 = op3.sites();  // Returns [2, 0]

Difference from support()

• sites() is only available for individual Spin-½ operators (e.g., X(0), XX(0,1), CNOT(2,0))
• support() works on individual operators, OperatorProduct, and OperatorSum
• sites() returns sites in definition order: CNOT(2,0).sites() returns [2, 0]
• support() returns an unordered set of all sites involved

Getting operator name

Use name() to get the operator name only as a string:

var op1 = X(0);
var name1 = name(op1);  // Returns 'X'

var op2 = XX(0, 1);
var name2 = name(op2);  // Returns 'XX'

var op3 = CNOT(0, 1);
var name3 = op3.name();  // Returns 'CNOT'

info

sites() and name() are designed to work only with individual Spin-½ operators (ID, X, Y, Z, XX, YY, CNOT, etc.).

Computing operator norm

By default, we use the spectral norm, i.e., the largest eigenvalue in absolute value of the operator. You can compute the norm of any operator using the norm() function:

var op1 = norm(X(0));           // Norm of a single Pauli operator
var op2 = norm(X(0) + Y(1));    // Norm of an `operator_sum`
var op3 = norm(2.0 * Z(0));     // Norm includes the coefficient

For a term with a coefficient and an operator, the norm behaves as you expect and combines the two:

$$|\text{coefficient}| \times \text{operator norm}$$

This effective norm is used by the prune transformations to determine which terms are negligible in an operator_sum.

For operator_sum we instead employ an effective norm that uses the triangle inequality,

$$|| \sum_m A_m || \leq \sum_m ||A_m||.$$

A number of known norms are encoded, and unitarity is exploited to ensure expensive calculations are avoided in these operations when possible.

Interacting with vectors

Operators can interact with appropriately sized vectors. We support * for out of place operations meaning we return a new vector and >> for in place operations, meaning we modify the vector in place without needing an auxiliary vector. To interact with spin $\frac{1}{2}$ operators your vectors should be a power of two with the appropriate number of degrees of freedom.

// 12 qubits, your support should be in range 0 to 11.  
var L = 12;
// random vector for the example
var vec = random([2**L])
// normalize it
normalize(vec)

var prod = RotY(pi/5,1)*RotX(pi/2,3)*Z(0)*X(5)*XX(5,11);
// modify vector in place
prod >> vec
// vec is kept the same, new vector constructed

var newvec = prod*vec;

tip

operator_prod and all local operators have efficient in place operations that don't require an expensive auxiliary, so >> operations will be the fastest, most efficient mode of operation.

Documentation Contributors

Jonathon Riddell

Eunji Yoo

Sebastien J Avakian