Skip to content

Add convenience function qml.math.ceil_log2 #8972

New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

Sign up for GitHub

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Merged
dwierichs merged 18 commits into from
Jan 25, 2026
Merged

Add convenience function qml.math.ceil_log2 #8972

Show file tree
Hide file tree
Changes from all commits
Commits
Show all changes
18 commits
Select commit Hold shift + click to select a range
dd14405
start
dwierichs Jan 22, 2026
2428f61
rank
dwierichs Jan 23, 2026
6ea98c1
use ceil_log2
dwierichs Jan 23, 2026
23ea7a6
reset binary stuff
dwierichs Jan 23, 2026
c34890f
mroe
dwierichs Jan 23, 2026
7c15f17
more
dwierichs Jan 23, 2026
2e9f81a
Merge branch 'master' into ceil-log2
dwierichs Jan 23, 2026
6d42e3e
last
dwierichs Jan 23, 2026
72d716f
debug
dwierichs Jan 23, 2026
8e11f7e
changelog
dwierichs Jan 23, 2026
8e61bf9
copy fix
dwierichs Jan 23, 2026
37e149f
copy fix
dwierichs Jan 23, 2026
b44ccb9
tiny
dwierichs Jan 23, 2026
f7a8ac4
salvage special cases using log2(0)=-inf
dwierichs Jan 23, 2026
34bc2a7
Merge branch 'master' into ceil-log2
dwierichs Jan 23, 2026
875578e
remaining occurrences
dwierichs Jan 23, 2026
bde49f7
remaining occurrences
dwierichs Jan 23, 2026
219bdd5
Merge branch 'master' into ceil-log2
dwierichs Jan 23, 2026
File filter

Filter by extension

Filter by extension
Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
5 changes: 5 additions & 0 deletions doc/releases/changelog-dev.md
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -2,6 +2,11 @@

<h3>New features since last release</h3>

* Added a convenience function :func:`~.math.ceil_log2` that computes the ceiling of the base-2
logarithm of its input and casts the result to an ``int``. It is equivalent to
``int(np.ceil(np.log2(n)))``.
[(#8972)](https://github.com/PennyLaneAI/pennylane/pull/8972)

* Added a ``qml.gate_sets`` that contains pre-defined gate sets such as ``qml.gate_sets.CLIFFORD_T_PLUS_RZ``
that can be plugged into the ``gate_set`` argument of the :func:`~pennylane.transforms.decompose` transform.
[(#8915)](https://github.com/PennyLaneAI/pennylane/pull/8915)
Expand Down
2 changes: 1 addition & 1 deletion pennylane/bose/bosonic_mapping.py
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -103,7 +103,7 @@ def _(bose_operator: BoseWord, n_states, tol=None):
raise ValueError(
f"Number of allowed bosonic states cannot be less than 2, provided {n_states}."
)
nqub_per_boson = int(np.ceil(np.log2(n_states)))
nqub_per_boson = math.ceil_log2(n_states)

creation = np.zeros((n_states, n_states))
for s in range(n_states - 1):
Expand Down
115 changes: 43 additions & 72 deletions pennylane/estimator/qpe_resources/first_quantization.py
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -19,6 +19,7 @@
import numpy as np
import scipy as sp

from pennylane.math import ceil_log2
from pennylane.operation import Operation


Expand Down Expand Up @@ -212,7 +213,7 @@ def success_prob(n, br):
if br <= 0 or not isinstance(br, int):
raise ValueError("br must be a positive integer.")

c = n / 2 ** np.ceil(np.log2(n))
c = n / 2 ** ceil_log2(n)
d = 2 * np.pi / 2**br

theta = d * np.round((1 / d) * np.arcsin(np.sqrt(1 / (4 * c))))
Expand Down Expand Up @@ -342,7 +343,7 @@ def norm(n, eta, omega, error, br=7, charge=0, cubic=True, vectors=None):
error_uv = 0.01 * error

# taken from Eq. (22) of PRX Quantum 2, 040332 (2021)
n_p = int(np.ceil(np.log2(n ** (1 / 3) + 1)))
n_p = ceil_log2(n ** (1 / 3) + 1)

n0 = n ** (1 / 3)
lambda_nu = ( # expression is taken from Eq. (F6) of PRX 8, 011044 (2018)
Expand Down Expand Up @@ -472,11 +473,11 @@ def unitary_cost(n, eta, omega, error, br=7, charge=0):
l_nu = 2 * np.pi * n ** (2 / 3)

# defined in the third and second paragraphs of page 15 of PRX Quantum 2, 040332 (2021)
n_eta = np.ceil(np.log2(eta))
n_etaz = np.ceil(np.log2(eta + 2 * l_z))
n_eta = ceil_log2(eta)
n_etaz = ceil_log2(eta + 2 * l_z)

# n_p is taken from Eq. (22)
n_p = int(np.ceil(np.log2(n ** (1 / 3) + 1)))
n_p = ceil_log2(n ** (1 / 3) + 1)

# errors in Eqs. (132-134) are set to be 0.01 of the algorithm error
error_t = alpha * error
Expand Down Expand Up @@ -667,7 +668,7 @@ def qubit_cost(n, eta, omega, error, br=7, charge=0, cubic=True, vectors=None):
l_nu = 2 * np.pi * n ** (2 / 3)

# n_p is taken from Eq. (22) of PRX Quantum 2, 040332 (2021)
n_p = np.ceil(np.log2(n ** (1 / 3) + 1))
n_p = ceil_log2(n ** (1 / 3) + 1)

# errors in Eqs. (132-134) of PRX Quantum 2, 040332 (2021),
# set to 0.01 of the algorithm error
Expand All @@ -693,8 +694,8 @@ def qubit_cost(n, eta, omega, error, br=7, charge=0, cubic=True, vectors=None):

# the expression for computing the cost is taken from Eq. (101) of arXiv:2204.11890v1
qubits = 3 * eta * n_p + 4 * n_m * n_p + 12 * n_p
qubits += 2 * (np.ceil(np.log2(np.ceil(np.pi * lamb / (2 * error_qpe))))) + 5 * n_m
qubits += 2 * np.ceil(np.log2(eta)) + 3 * n_p**2 + np.ceil(np.log2(eta + 2 * l_z))
qubits += 2 * (ceil_log2(np.ceil(np.pi * lamb / (2 * error_qpe)))) + 5 * n_m
qubits += 2 * ceil_log2(eta) + 3 * n_p**2 + ceil_log2(eta + 2 * l_z)
qubits += np.maximum(5 * n_p + 1, 5 * n_r - 4) + np.maximum(n_t, n_r + 1) + 33

return int(np.ceil(qubits))
Expand Down Expand Up @@ -746,23 +747,19 @@ def _norm_noncubic(
error_uv = alpha * error

# taken from Eq. (22) of PRX Quantum 2, 040332 (2021)
n_p = int(np.ceil(np.log2(n ** (1 / 3) + 1)))
n_p = ceil_log2(n ** (1 / 3) + 1)

n0 = n ** (1 / 3)

# defined in Eq. (F3) of arXiv:2302.07981v1 (2023)
bmin = np.min(np.linalg.svd(recip_vectors)[1])

n_m = int(
np.ceil(
np.log2( # taken from Eq. (132) of PRX Quantum 2, 040332 (2021) with
# modifications taken from arXiv:2302.07981v1 (2023)
(8 * np.pi * eta)
/ (error_uv * omega * bmin**2)
* (eta - 1 + 2 * l_z)
* (7 * 2 ** (n_p + 1) - 9 * n_p - 11 - 3 * 2 ** (-1 * n_p))
)
)
n_m = ceil_log2( # taken from Eq. (132) of PRX Quantum 2, 040332 (2021) with
# modifications taken from arXiv:2302.07981v1 (2023)
(8 * np.pi * eta)
/ (error_uv * omega * bmin**2)
* (eta - 1 + 2 * l_z)
* (7 * 2 ** (n_p + 1) - 9 * n_p - 11 - 3 * 2 ** (-1 * n_p))
)

lambda_nu = ( # expression is taken from Eq. (F6) of PRX 8, 011044 (2018)
Expand Down Expand Up @@ -863,7 +860,7 @@ def _qubit_cost_noncubic(n, eta, error, br, charge, vectors):
l_nu = 2 * np.pi * n ** (2 / 3)

# taken from Eq. (22) of PRX Quantum 2, 040332 (2021)
n_p = np.ceil(np.log2(n ** (1 / 3) + 1))
n_p = ceil_log2(n ** (1 / 3) + 1)

# defined in Eq. (F3) of arXiv:2302.07981v1 (2023)
bmin = np.min(np.linalg.svd(recip_vectors)[1])
Expand All @@ -874,17 +871,13 @@ def _qubit_cost_noncubic(n, eta, error, br, charge, vectors):
error_t, error_r, error_m, error_b = [alpha * error] * 4

# parameters taken from PRX Quantum 2, 040332 (2021)
n_t = int(np.ceil(np.log2(np.pi * lambda_total / error_t))) # Eq. (134)
n_r = int(np.ceil(np.log2((eta * l_z * l_nu) / (error_r * omega ** (1 / 3))))) # Eq. (133)
n_m = int( # taken from Eq. (J13) arXiv:2302.07981v1 (2023)
np.ceil(
np.log2(
(8 * np.pi * eta)
/ (error_m * omega * bmin**2)
* (eta - 1 + 2 * l_z)
* (7 * 2 ** (n_p + 1) - 9 * n_p - 11 - 3 * 2 ** (-1 * n_p))
)
)
n_t = ceil_log2(np.pi * lambda_total / error_t) # Eq. (134)
n_r = ceil_log2((eta * l_z * l_nu) / (error_r * omega ** (1 / 3))) # Eq. (133)
n_m = ceil_log2( # taken from Eq. (J13) arXiv:2302.07981v1 (2023)
(8 * np.pi * eta)
/ (error_m * omega * bmin**2)
* (eta - 1 + 2 * l_z)
* (7 * 2 ** (n_p + 1) - 9 * n_p - 11 - 3 * 2 ** (-1 * n_p))
)

# qpe error obtained to satisfy a modification of
Expand All @@ -894,33 +887,23 @@ def _qubit_cost_noncubic(n, eta, error, br, charge, vectors):

# adapted from equations (L1, L2) in Appendix L of arXiv:2302.07981v1 (2023)
clean_temp_H_cost = max([5 * n_r - 4, 5 * n_p + 1]) + max([5, n_m + 3 * n_p])
reflection_cost = (
np.ceil(np.log2(eta + 2 * l_z)) + 2 * np.ceil(np.log2(eta)) + 6 * n_p + n_m + 16 + 3
)
reflection_cost = ceil_log2(eta + 2 * l_z) + 2 * ceil_log2(eta) + 6 * n_p + n_m + 16 + 3
clean_temp_cost = max([clean_temp_H_cost, reflection_cost])

# the expression for computing the cost is taken from Appendix C of
# PRX Quantum 2, 040332 (2021) and adapting to non-cubic using Appendix L of
# arXiv:2302.07981v1 (2023)
clean_cost = 3 * eta * n_p
clean_cost += np.ceil(np.log2(np.ceil(np.pi * lambda_total / (2 * error_qpe))))
clean_cost += 1 + 1 + np.ceil(np.log2(eta + 2 * l_z)) + 3 + 3
clean_cost += 2 * np.ceil(np.log2(eta)) + 5 + 3 * (n_p + 1)
clean_cost += ceil_log2(np.ceil(np.pi * lambda_total / (2 * error_qpe)))
clean_cost += 1 + 1 + ceil_log2(eta + 2 * l_z) + 3 + 3
clean_cost += 2 * ceil_log2(eta) + 5 + 3 * (n_p + 1)
clean_cost += n_p + n_m + 3 * n_p + 2 + 2 * n_p + 1 + 1 + 2 + 2 * n_p + 6 + 1

clean_cost += clean_temp_cost

# taken from Eq. (J7) of arXiv:2302.07981v1 (2023)
n_b = np.ceil(
np.log2(
4
* np.pi
* eta
* 2 ** (2 * n_p - 2)
* np.abs(recip_vectors @ recip_vectors.T).flatten().sum()
/ error_b
)
)
sum_abs = np.abs(recip_vectors @ recip_vectors.T).flatten().sum()
n_b = ceil_log2(4 * np.pi * eta * 2 ** (2 * n_p - 2) * sum_abs / error_b)
clean_cost += np.max([n_r + 1, n_t, n_b]) + 6 + n_m + 1

return int(np.ceil(clean_cost))
Expand Down Expand Up @@ -967,32 +950,28 @@ def _unitary_cost_noncubic(n, eta, error, br, charge, vectors):
l_nu = 2 * np.pi * n ** (2 / 3)

# defined in the third and second paragraphs of page 15 of PRX Quantum 2, 040332 (2021)
n_eta = np.ceil(np.log2(eta))
n_etaz = np.ceil(np.log2(eta + 2 * l_z))
n_eta = ceil_log2(eta)
n_etaz = ceil_log2(eta + 2 * l_z)

# taken from Eq. (22) of PRX Quantum 2, 040332 (2021)
n_p = int(np.ceil(np.log2(n ** (1 / 3) + 1)))
n_p = ceil_log2(n ** (1 / 3) + 1)

# errors in Eqs. (132-134) of PRX Quantum 2, 040332 (2021)
error_t, error_r, error_m = [alpha * error] * 3

# parameters taken from PRX Quantum 2, 040332 (2021)
n_t = int(np.ceil(np.log2(np.pi * lambda_total / error_t))) # Eq. (134)
n_r = int(np.ceil(np.log2((eta * l_z * l_nu) / (error_r * omega ** (1 / 3))))) # Eq. (133)
n_t = ceil_log2(np.pi * lambda_total / error_t) # Eq. (134)
n_r = ceil_log2((eta * l_z * l_nu) / (error_r * omega ** (1 / 3))) # Eq. (133)

# defined in Eq. (F3) of arXiv:2302.07981v1 (2023)
bmin = np.min(np.linalg.svd(recip_vectors)[1])

# equivalent to Eq. (J13) of arXiv:2302.07981v1 (2023)
n_m = int(
np.ceil(
np.log2(
(8 * np.pi * eta)
/ (error_m * omega * bmin**2)
* (eta - 1 + 2 * l_z)
* (7 * 2 ** (n_p + 1) - 9 * n_p - 11 - 3 * 2 ** (-1 * n_p))
)
)
n_m = ceil_log2(
(8 * np.pi * eta)
/ (error_m * omega * bmin**2)
* (eta - 1 + 2 * l_z)
* (7 * 2 ** (n_p + 1) - 9 * n_p - 11 - 3 * 2 ** (-1 * n_p))
)

e_r = FirstQuantization._cost_qrom(l_z)
Expand All @@ -1006,16 +985,8 @@ def _unitary_cost_noncubic(n, eta, error, br, charge, vectors):

# taken from Eq. (J7) of arXiv:2302.07981v1 (2023)
error_b = alpha * error
n_b = np.ceil(
np.log2(
2
* np.pi
* eta
* 2 ** (2 * n_p - 2)
* np.abs(recip_vectors @ recip_vectors.T).flatten().sum()
/ error_b
)
)
sum_abs = np.abs(recip_vectors @ recip_vectors.T).flatten().sum()
n_b = ceil_log2(2 * np.pi * eta * 2 ** (2 * n_p - 2) * sum_abs / error_b)

n_dirty = FirstQuantization._qubit_cost_noncubic(n, eta, error, br, charge, vectors)
ms_cost = FirstQuantization._momentum_state_qrom(n_p, n_m, n_dirty, n_tof, kappa=1)[0]
Expand Down Expand Up @@ -1048,7 +1019,7 @@ def _momentum_state_qrom(n_p, n_m, n_dirty, n_tof, kappa):
else:
beta_gate = max([np.floor(2 * x / (3 * n_m / kappa) * np.log(2)), 1])
beta = np.min([beta_dirty, beta_gate, beta_parallel])
ms_cost_qrom = 2 * np.ceil(x / beta) + 3 * np.ceil(n_m / kappa) * np.ceil(np.log2(beta))
ms_cost_qrom = 2 * np.ceil(x / beta) + 3 * np.ceil(n_m / kappa) * ceil_log2(beta)

ms_cost = 2 * ms_cost_qrom + n_m + 8 * (n_p - 1) + 6 * n_p + 2 + 2 * n_p + n_m + 2

Expand Down
15 changes: 8 additions & 7 deletions pennylane/estimator/qpe_resources/second_quantization.py
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -18,6 +18,7 @@
# pylint: disable=no-self-use, too-many-arguments, too-many-instance-attributes, too-many-positional-arguments
import numpy as np

from pennylane.math import ceil_log2
from pennylane.operation import Operation
from pennylane.qchem import factorize

Expand Down Expand Up @@ -372,9 +373,9 @@ def unitary_cost(n, rank_r, rank_m, rank_max, br=7, alpha=10, beta=20):
eta = np.array([np.log2(n) for n in range(1, rank_r + 1) if rank_r % n == 0])
eta = int(np.max([n for n in eta if n % 1 == 0]))

nxi = np.ceil(np.log2(rank_max)) # Eq. (C14) of PRX Quantum 2, 030305 (2021)
nl = np.ceil(np.log2(rank_r + 1)) # Eq. (C14) of PRX Quantum 2, 030305 (2021)
nlxi = np.ceil(np.log2(rank_rm + n / 2)) # Eq. (C15) of PRX Quantum 2, 030305 (2021)
nxi = ceil_log2(rank_max) # Eq. (C14) of PRX Quantum 2, 030305 (2021)
nl = ceil_log2(rank_r + 1) # Eq. (C14) of PRX Quantum 2, 030305 (2021)
nlxi = ceil_log2(rank_rm + n / 2) # Eq. (C15) of PRX Quantum 2, 030305 (2021)

bp1 = nl + alpha # Eq. (C27) of PRX Quantum 2, 030305 (2021)
bo = nxi + nlxi + br + 1 # Eq. (C29) of PRX Quantum 2, 030305 (2021)
Expand Down Expand Up @@ -531,9 +532,9 @@ def qubit_cost(n, lamb, error, rank_r, rank_m, rank_max, br=7, alpha=10, beta=20

rank_rm = rank_r * rank_m

nxi = np.ceil(np.log2(rank_max)) # Eq. (C14) of PRX Quantum 2, 030305 (2021)
nl = np.ceil(np.log2(rank_r + 1)) # Eq. (C14) of PRX Quantum 2, 030305 (2021)
nlxi = np.ceil(np.log2(rank_rm + n / 2)) # Eq. (C15) of PRX Quantum 2, 030305 (2021)
nxi = ceil_log2(rank_max) # Eq. (C14) of PRX Quantum 2, 030305 (2021)
nl = ceil_log2(rank_r + 1) # Eq. (C14) of PRX Quantum 2, 030305 (2021)
nlxi = ceil_log2(rank_rm + n / 2) # Eq. (C15) of PRX Quantum 2, 030305 (2021)

bo = nxi + nlxi + br + 1 # Eq. (C29) of PRX Quantum 2, 030305 (2021)
bp2 = nxi + alpha + 2 # Eq. (C31) of PRX Quantum 2, 030305 (2021)
Expand All @@ -543,7 +544,7 @@ def qubit_cost(n, lamb, error, rank_r, rank_m, rank_max, br=7, alpha=10, beta=20
# the cost is computed using Eq. (C40) of PRX Quantum 2, 030305 (2021)
e_cost = DoubleFactorization.estimation_cost(lamb, error)
cost = n + 2 * nl + nxi + 3 * alpha + beta + bo + bp2
cost += kr * n * beta / 2 + 2 * np.ceil(np.log2(e_cost + 1)) + 7
cost += kr * n * beta / 2 + 2 * ceil_log2(e_cost + 1) + 7

return int(cost)

Expand Down
28 changes: 8 additions & 20 deletions pennylane/estimator/templates/qubitize.py
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -12,10 +12,6 @@
# See the License for the specific language governing permissions and
# limitations under the License.
r"""Resource operators for PennyLane subroutine templates."""
import math

import numpy as np

from pennylane.estimator.compact_hamiltonian import THCHamiltonian
from pennylane.estimator.ops.op_math.controlled_ops import MultiControlledX, Toffoli
from pennylane.estimator.ops.op_math.symbolic import Adjoint, Controlled
Expand All @@ -30,6 +26,7 @@
from pennylane.estimator.templates.select import SelectTHC
from pennylane.estimator.templates.stateprep import PrepTHC
from pennylane.estimator.wires_manager import Allocate, Deallocate
from pennylane.math import ceil_log2
from pennylane.wires import Wires, WiresLike

# pylint: disable=signature-differs, arguments-differ, too-many-arguments
Expand Down Expand Up @@ -125,7 +122,7 @@ def __init__(
num_orb = thc_ham.num_orbitals
tensor_rank = thc_ham.tensor_rank
num_coeff = num_orb + tensor_rank * (tensor_rank + 1) / 2 # N+M(M+1)/2
coeff_register = int(math.ceil(math.log2(num_coeff)))
coeff_register = ceil_log2(num_coeff)

if coeff_precision is None:
coeff_precision = prep_op.coeff_precision if prep_op else 15
Expand Down Expand Up @@ -156,13 +153,8 @@ def __init__(
# The total algorithmic qubits are thus given by: N + 2*n_M + ceil(log(d)) + \aleph + 6 + m
# where \aleph is coeff_precision, m = 2n_M + \aleph + 2, N = 2*num_orb,
# d = num_orb + tensor_rank(tensor_rank+1)/2, and n_M = log_2(tensor_rank+1).
self.num_wires = (
num_orb * 2
+ 4 * int(np.ceil(math.log2(tensor_rank + 1)))
+ coeff_register
+ 8
+ coeff_precision
)
n_M = ceil_log2(tensor_rank + 1)
self.num_wires = num_orb * 2 + 4 * n_M + coeff_register + 8 + coeff_precision
if wires is not None and len(Wires(wires)) != self.num_wires:
raise ValueError(f"Expected {self.num_wires} wires, got {len(Wires(wires))}")
super().__init__(wires=wires)
Expand Down Expand Up @@ -232,19 +224,15 @@ def resource_rep(
num_orb = thc_ham.num_orbitals
tensor_rank = thc_ham.tensor_rank
num_coeff = num_orb + tensor_rank * (tensor_rank + 1) / 2 # N+M(M+1)/2
coeff_register = int(math.ceil(math.log2(num_coeff)))
coeff_register = ceil_log2(num_coeff)

if coeff_precision is None:
coeff_precision = prep_op.params["coeff_precision"] if prep_op else 15
if rotation_precision is None:
rotation_precision = select_op.params["rotation_precision"] if select_op else 15

num_wires = (
num_orb * 2
+ 4 * int(np.ceil(math.log2(tensor_rank + 1)))
+ coeff_register
+ 8
+ coeff_precision
num_orb * 2 + 4 * ceil_log2(tensor_rank + 1) + coeff_register + 8 + coeff_precision
)

params = {
Expand Down Expand Up @@ -297,7 +285,7 @@ def resource_decomp(
gate_list = []

tensor_rank = thc_ham.tensor_rank
m_register = int(np.ceil(np.log2(tensor_rank)))
m_register = ceil_log2(tensor_rank)

select_kwargs = {
"thc_ham": thc_ham,
Expand Down Expand Up @@ -367,7 +355,7 @@ def controlled_resource_decomp(
rotation_precision = target_resource_params.get("rotation_precision")

tensor_rank = thc_ham.tensor_rank
m_register = int(np.ceil(np.log2(tensor_rank)))
m_register = ceil_log2(tensor_rank)

if num_ctrl_wires > 1:
mcx = resource_rep(
Expand Down
Loading