qml.resource.DoubleFactorization

class DoubleFactorization(one_electron, two_electron, error=0.0016, rank_r=None, rank_m=None, rank_max=None, tol_factor=1e-05, tol_eigval=1e-05, br=7, alpha=10, beta=20, chemist_notation=False)[source]

Bases: pennylane.operation.Operation

Estimate the number of non-Clifford gates and logical qubits for a quantum phase estimation algorithm in second quantization with a double-factorized Hamiltonian.

Atomic units are used throughout the class.

Parameters
  • one_electron (array[array[float]]) – one-electron integrals

  • two_electron (tensor_like) – two-electron integrals

  • error (float) – target error in the algorithm

  • rank_r (int) – rank of the first factorization of the two-electron integral tensor

  • rank_m (int) – average rank of the second factorization of the two-electron integral tensor

  • tol_factor (float) – threshold error value for discarding the negligible factors

  • tol_eigval (float) – threshold error value for discarding the negligible factor eigenvalues

  • br (int) – number of bits for ancilla qubit rotation

  • alpha (int) – number of bits for the keep register

  • beta (int) – number of bits for the rotation angles

  • chemist_notation (bool) – if True, the two-electron integrals need to be in chemist notation

Example

>>> symbols  = ['O', 'H', 'H']
>>> geometry = np.array([[0.00000000,  0.00000000,  0.28377432],
>>>                      [0.00000000,  1.45278171, -1.00662237],
>>>                      [0.00000000, -1.45278171, -1.00662237]], requires_grad = False)
>>> mol = qml.qchem.Molecule(symbols, geometry, basis_name='sto-3g')
>>> core, one, two = qml.qchem.electron_integrals(mol)()
>>> algo = DoubleFactorization(one, two)
>>> print(algo.lamb,  # the 1-Norm of the Hamiltonian
>>>       algo.gates, # estimated number of non-Clifford gates
>>>       algo.qubits # estimated number of logical qubits
>>>       )
53.62085493277858 103969925 290

To estimate the gate and qubit costs for implementing this method, the Hamiltonian needs to be factorized using the factorize() function following [PRX Quantum 2, 030305 (2021)]. The objective of the factorization is to find a set of symmetric matrices, \(L^{(r)}\), such that the two-electron integral tensor in chemist notation, \(V\), can be computed as

\[V_{ijkl} = \sum_r^R L_{ij}^{(r)} L_{kl}^{(r) T},\]

with the rank \(R \leq n^2\), where \(n\) is the number of molecular orbitals. The matrices \(L^{(r)}\) are diagonalized and for each matrix the eigenvalues that are smaller than a given threshold (and their corresponding eigenvectors) are discarded. The average number of the retained eigenvalues, \(M\), determines the rank of the second factorization step. The 1-norm of the Hamiltonian can then be computed using the norm() function from the electron integrals and the eigenvalues of the matrices \(L^{(r)}\).

The total number of gates and qubits for implementing the quantum phase estimation algorithm for the given Hamiltonian can then be computed using the functions gate_cost() and qubit_cost() with a target error that has the default value of 0.0016 Ha (chemical accuracy). The costs are computed using Eqs. (C39-C40) of [PRX Quantum 2, 030305 (2021)].

grad_method

num_wires

grad_method = None
num_wires = -1

estimation_cost(lamb, error)

Return the number of calls to the unitary needed to achieve the desired error in quantum phase estimation.

gate_cost(n, lamb, error, rank_r, rank_m, …)

Return the total number of Toffoli gates needed to implement the double factorization algorithm.

norm(one, two, eigvals)

Return the 1-norm of a molecular Hamiltonian from the one- and two-electron integrals and eigenvalues of the factorized two-electron integral tensor.

qubit_cost(n, lamb, error, rank_r, rank_m, …)

Return the number of logical qubits needed to implement the double factorization method.

unitary_cost(n, rank_r, rank_m, rank_max[, …])

Return the number of Toffoli gates needed to implement the qubitization unitary operator for the double factorization algorithm.

static estimation_cost(lamb, error)[source]

Return the number of calls to the unitary needed to achieve the desired error in quantum phase estimation.

The expression for computing the cost is taken from Eq. (45) of [PRX Quantum 2, 030305 (2021)].

Parameters
  • lamb (float) – 1-norm of a second-quantized Hamiltonian

  • error (float) – target error in the algorithm

Returns

number of calls to unitary

Return type

int

Example

>>> lamb = 72.49779513025341
>>> error = 0.001
>>> estimation_cost(lamb, error)
113880
static gate_cost(n, lamb, error, rank_r, rank_m, rank_max, br=7, alpha=10, beta=20)[source]

Return the total number of Toffoli gates needed to implement the double factorization algorithm.

The expression for computing the cost is taken from Eqs. (45) and (C39) of [PRX Quantum 2, 030305 (2021)].

Parameters
  • n (int) – number of molecular spin-orbitals

  • lamb (float) – 1-norm of a second-quantized Hamiltonian

  • error (float) – target error in the algorithm

  • rank_r (int) – rank of the first factorization of the two-electron integral tensor

  • rank_m (float) – average rank of the second factorization of the two-electron tensor

  • rank_max (int) – maximum rank of the second factorization of the two-electron tensor

  • br (int) – number of bits for ancilla qubit rotation

  • alpha (int) – number of bits for the keep register

  • beta (int) – number of bits for the rotation angles

Returns

the number of Toffoli gates for the double factorization method

Return type

int

Example

>>> n = 14
>>> lamb = 52.98761457453095
>>> error = 0.001
>>> rank_r = 26
>>> rank_m = 5.5
>>> rank_max = 7
>>> br = 7
>>> alpha = 10
>>> beta = 20
>>> gate_cost(n, lamb, error, rank_r, rank_m, rank_max, br, alpha, beta)
167048631
static norm(one, two, eigvals)[source]

Return the 1-norm of a molecular Hamiltonian from the one- and two-electron integrals and eigenvalues of the factorized two-electron integral tensor.

The 1-norm of a double-factorized molecular Hamiltonian is computed using Eqs. (15-17) of [Phys. Rev. Research 3, 033055 (2021)]

\[\lambda = ||T|| + \frac{1}{4} \sum_r ||L^{(r)}||^2,\]

where the Schatten 1-norm for a given matrix \(T\) is defined as

\[||T|| = \sum_k |\text{eigvals}[T]_k|.\]

The matrices \(L^{(r)}\) are obtained from factorization of the two-electron integral tensor \(V\) such that

\[V_{ijkl} = \sum_r L_{ij}^{(r)} L_{kl}^{(r) T}.\]

The matrix \(T\) is constructed from the one- and two-electron integrals as

\[T = h_{ij} - \frac{1}{2} \sum_l V_{illj} + \sum_l V_{llij}.\]

Note that the two-electron integral tensor must be arranged in chemist notation.

Parameters
  • one (array[array[float]]) – one-electron integrals

  • two (array[array[float]]) – two-electron integrals

  • eigvals (array[float]) – eigenvalues of the matrices obtained from factorizing the two-electron integral tensor

Returns

1-norm of the Hamiltonian

Return type

array[float]

Example

>>> symbols  = ['H', 'H', 'O']
>>> geometry = np.array([[0.00000000,  0.00000000,  0.28377432],
>>>                      [0.00000000,  1.45278171, -1.00662237],
>>>                      [0.00000000, -1.45278171, -1.00662237]], requires_grad=False)
>>> mol = qml.qchem.Molecule(symbols, geometry, basis_name='sto-3g')
>>> core, one, two = qml.qchem.electron_integrals(mol)()
>>> two = np.swapaxes(two, 1, 3) # convert to the chemists notation
>>> _, eigvals, _ = qml.qchem.factorize(two, 1e-5)
>>> print(norm(one, two, eigvals))
52.98762043980203
static qubit_cost(n, lamb, error, rank_r, rank_m, rank_max, br=7, alpha=10, beta=20)[source]

Return the number of logical qubits needed to implement the double factorization method.

The expression for computing the cost is taken from Eq. (C40) of [PRX Quantum 2, 030305 (2021)].

Parameters
  • n (int) – number of molecular spin-orbitals

  • lamb (float) – 1-norm of a second-quantized Hamiltonian

  • error (float) – target error in the algorithm

  • rank_r (int) – rank of the first factorization of the two-electron integral tensor

  • rank_m (float) – average rank of the second factorization of the two-electron tensor

  • rank_max (int) – maximum rank of the second factorization of the two-electron tensor

  • br (int) – number of bits for ancilla qubit rotation

  • alpha (int) – number of bits for the keep register

  • beta (int) – number of bits for the rotation angles

Returns

number of logical qubits for the double factorization method

Return type

int

Example

>>> n = 14
>>> lamb = 52.98761457453095
>>> error = 0.001
>>> rank_r = 26
>>> rank_m = 5.5
>>> rank_max = 7
>>> br = 7
>>> alpha = 10
>>> beta = 20
>>> qubit_cost(n, lamb, error, rank_r, rank_m, rank_max, br, alpha, beta)
292
static unitary_cost(n, rank_r, rank_m, rank_max, br=7, alpha=10, beta=20)[source]

Return the number of Toffoli gates needed to implement the qubitization unitary operator for the double factorization algorithm.

The expression for computing the cost is taken from Eq. (C39) of [PRX Quantum 2, 030305 (2021)].

Parameters
  • n (int) – number of molecular spin-orbitals

  • rank_r (int) – rank of the first factorization of the two-electron integral tensor

  • rank_m (float) – average rank of the second factorization of the two-electron tensor

  • rank_max (int) – maximum rank of the second factorization of the two-electron tensor

  • br (int) – number of bits for ancilla qubit rotation

  • alpha (int) – number of bits for the keep register

  • beta (int) – number of bits for the rotation angles

Returns

number of Toffoli gates to implement the qubitization unitary

Return type

int

Example

>>> n = 14
>>> rank_r = 26
>>> rank_m = 5.5
>>> rank_max = 7
>>> br = 7
>>> alpha = 10
>>> beta = 20
>>> unitary_cost(n, rank_r, rank_m, rank_max, br, alpha, beta)
2007