Reduced Density Matrices
Reduced Density Matrices (RDMs) are a use case where measurement reduction in
combination with computables is very useful. Naively, we can evaluate each term in the n-body RDM
by evaluating many ExpectationValue
computables independently
on each spin-traced excitation operator.
However, across the entire RDM, this would result in many redundant
measurements, making it advantageous to create a single, composite computable for an RDM ndarray
and use a protocol to build the optimal number of circuits for the whole RDM.
To build a 1-RDM operator, we can use InQuanto’s FermionSpace
:
import numpy
from inquanto.spaces import FermionSpace
space = FermionSpace(4)
rdm_operators = space.construct_one_body_spatial_rdm_operators((FermionSpace.SPIN_UP, FermionSpace.SPIN_UP))
print(rdm_operators)
print(type(rdm_operators))
[[{((0, 1), (0, 0)): 1.0} {((0, 1), (2, 0)): 1.0}]
[{((2, 1), (0, 0)): 1.0} {((2, 1), (2, 0)): 1.0}]]
<class 'numpy.ndarray'>
As a result, we obtain a numpy ndarray
of FermionOperator
objects.
To compute the expectation values, we also need a specific state.
In this example, we will use the TrotterAnsatz
:
from inquanto.ansatzes import TrotterAnsatz
from inquanto.operators import QubitOperator, QubitOperatorList
from inquanto.states import QubitState
ansatz = TrotterAnsatz(
QubitOperatorList.from_list([QubitOperator("Y0 X1 X2 X3", 1j)]),
QubitState([1, 1, 0, 0]),
)
parameters = ansatz.state_symbols.construct_from_array([0.11])
With the ansatz and RDM operators, we can construct a single computable for the RDM:
from inquanto.computables import ExpectationValueNonHermitian
from inquanto.computables.primitive import ComputableNDArray
# We make a broadcast function over the ndarray that converts every element into
# an ExpectationValueNonHermitian computable
ndarray_broadcast = numpy.vectorize(
lambda operator: ExpectationValueNonHermitian(ansatz, operator.qubit_encode())
)
# Then, we can convert the ndarray computables to a ComputableNDArray of computables
# to ensure that we can run and evaluate recursively.
rdm_computable = ComputableNDArray(ndarray_broadcast(rdm_operators))
print(rdm_computable)
print(type(rdm_computable))
[[ExpectationValueNonHermitian(state=<inquanto.ansatzes._trotter_ansatz.TrotterAnsatz, qubits=4, gates=33, symbols=1>, kernel={(): 0.5, (Zq[0]): -0.5})
ExpectationValueNonHermitian(state=<inquanto.ansatzes._trotter_ansatz.TrotterAnsatz, qubits=4, gates=33, symbols=1>, kernel={(Yq[0], Zq[1], Xq[2]): (-0-0.25j), (Yq[0], Zq[1], Yq[2]): 0.25, (Xq[0], Zq[1], Xq[2]): 0.25, (Xq[0], Zq[1], Yq[2]): 0.25j})]
[ExpectationValueNonHermitian(state=<inquanto.ansatzes._trotter_ansatz.TrotterAnsatz, qubits=4, gates=33, symbols=1>, kernel={(Yq[0], Zq[1], Xq[2]): 0.25j, (Xq[0], Zq[1], Xq[2]): 0.25, (Yq[0], Zq[1], Yq[2]): 0.25, (Xq[0], Zq[1], Yq[2]): (-0-0.25j)})
ExpectationValueNonHermitian(state=<inquanto.ansatzes._trotter_ansatz.TrotterAnsatz, qubits=4, gates=33, symbols=1>, kernel={(): 0.5, (Zq[2]): -0.5})]]
<class 'inquanto.computables.primitive.collections._ndarray.ComputableNDArray'>
With this RDM computable, we can use protocols to create an evaluator function. We first use
default_evaluate()
to perform a
statevector calculation as a test:
print(rdm_computable.default_evaluate(parameters))
[[(0.9879487246653036+0j) (3.448477954431592e-18-1.3465787677460238e-19j)]
[(3.448477954431592e-18+1.3465787677460238e-19j)
(0.012051275334697197+0j)]]
Now, we use the shot-based PauliAveraging
protocol to
evaluate the RDM computable:
from pytket.extensions.qiskit import AerBackend
from inquanto.protocols import PauliAveraging
from pytket.partition import PauliPartitionStrat
# We can use the PauliAveraging protocol to get an evaluator
evaluator = (
PauliAveraging(
AerBackend(),
shots_per_circuit = 10000,
pauli_partition_strategy=PauliPartitionStrat.CommutingSets
)
.build_from(parameters, rdm_computable)
.run(seed=0)
.get_evaluator()
)
rdm_measured = rdm_computable.evaluate(evaluator)
print(rdm_measured)
[[(0.9867+0j) (0.0029+0.0041j)]
[(0.0029-0.0041j) (0.013299999999999979+0j)]]
InQuanto provides built-in support for RDM computables such as the above, for example with the
RestrictedOneBodyRDMComputable
class, which evaluates to a
RestrictedOneBodyRDM
. A suite of other density matrix composite
computables are supported,
for instance, the SpinlessNBodyRDMArrayRealComputable
can build n-body RDMs and reduce measurement circuits further by taking symmetries
into account:
from inquanto.ansatzes import FermionSpaceAnsatzChemicallyAwareUCCSD
from inquanto.computables.composite import SpinlessNBodyRDMArrayRealComputable
from inquanto.express import load_h5
from inquanto.mappings import QubitMappingJordanWigner
from inquanto.spaces import QubitSpace, FermionSpace
from inquanto.states import FermionState
qubit_hamiltonian = load_h5(
"h2_631g_symmetry.h5", as_tuple=True
).hamiltonian_operator.qubit_encode()
space = FermionSpace(
8,
point_group="D2h",
orb_irreps=numpy.asarray(["Ag", "Ag", "B1u", "B1u", "Ag", "Ag", "B1u", "B1u"]),
)
ansatz = FermionSpaceAnsatzChemicallyAwareUCCSD(
fermion_space=space, fermion_state=FermionState([1, 1, 0, 0, 0, 0, 0, 0])
)
qubit_space = QubitSpace(space.n_spin_orb)
symmetry_operators = qubit_space.symmetry_operators_z2(qubit_hamiltonian)
rdm1 = SpinlessNBodyRDMArrayRealComputable(
1,
space,
ansatz,
QubitMappingJordanWigner(),
symmetry_operators,
)
print(rdm1)
<inquanto.computables.composite.rdm._rdm_real.SpinlessNBodyRDMArrayRealComputable object at 0x7f6316185c90>
Above, we used the express module to load a Hamiltonian, and use a Fermion space
to generate the Z2 symmetries of our Hamiltonian.
The SpinlessNBodyRDMArrayRealComputable
uses the symmetry information to filter for Z2 symmetry.