discretelog

name: discretelog description: "Use when users ask about solving the discrete logarithm problem g^x ≡ y (mod P) with Shor's quantum algorithm, building/explaining DLP circuits, running simulator demos, or debugging post-processing (continued fractions, order recovery, congruence solving). Triggers: discrete log, DLP, Shor discrete logarithm, g^x mod P, modular exponentiation, continued fractions, quantum cryptography demo."

Discrete Logarithm Algorithm (DLG)

Purpose

Solves the discrete logarithm problem (DLP): given $g$, $y$, and prime $P$ with $g^x \equiv y \pmod{P}$, find $x$. The quantum algorithm runs in polynomial time $O(n^3)$ where $n = \log_2 P$, compared to the best classical sub-exponential algorithms.

Use this skill when you need to:

• Demonstrate a quantum attack on DLP-based cryptography (Diffie-Hellman, ECC).
• Understand the 2D quantum period-finding extension of Shor's algorithm.

Overview

1. Prepare three registers: reg_A ($n_{\text{count}}$ bits), reg_B ($n_{\text{count}}$ bits), work ($n_{\text{work}}$ bits).
2. Apply Hadamard to both A and B; initialize work to $|1\rangle$.
3. Apply controlled $g^a \bmod P$ (reg_A controls $\times g^{2^i}$) to work.
4. Apply controlled $(y^{-1})^b \bmod P$ (reg_B controls $\times (y^{-1})^{2^j}$) to work.
5. Apply IQFT to both A and B.
6. Measure A and B to get integers $(u, v)$.
7. Classical post-processing: continued fractions on $u/N$ to extract period $r$ and random $s$, then solve $sx \equiv \text{round}(v \cdot r/N) \pmod{r}$ for $x$.

Prerequisites

• Quantum Fourier Transform and QPE.
• Modular arithmetic, multiplicative order, modular inverse.
• Continued fractions algorithm.
• Python: numpy, Circuit, Register, unitarylab.library.IQFT.

Using the Provided Implementation

from unitarylab_algorithms import DiscreteLogAlgorithm

algo = DiscreteLogAlgorithm()
result = algo.run(
    g=3,   # Base
    y=6,   # Target: 3^x ≡ 6 (mod 7)
    P=7,   # Modulus (prime)
    backend='torch'
)

print(result['Found x'])                  # Found discrete log x
print(result['status'])                   # 'ok' on success, 'failed' otherwise
print(result['circuit_path'])             # SVG circuit diagram path
print(result.get('plot', []))             # List of output file dicts [{"format": ..., "filename": ...}]
print(result.get('Detected period r'))    # Detected group order r
print(result.get('Computation time (s)')) # Simulation time in seconds

Core Parameters Explained

Constructor DiscreteLogAlgorithm(text_mode, algo_dir):

run(g, y, P, backend, device, dtype) parameters:

Common misunderstandings:

• Both g and y must be coprime to P. This is validated; a ValueError is raised otherwise.
• The circuit uses $2 \cdot n_{\text{count}} + n_{\text{work}} = 2 \cdot 2\lfloor\log_2 P\rfloor + \lfloor\log_2 P\rfloor$ total qubits.
• Success is probabilistic; the algorithm relies on the measurement giving a valid continued fraction.

Return Fields

Implementation Architecture

DiscreteLogAlgorithm in algorithm.py structures the DLP algorithm into five ordered stages within run(), using a matrix-based modular multiplication oracle and a multi-step classical post-processing routine.

run(g, y, P, backend, device, dtype) — Five Stages:

Helper Methods:

• _get_modular_matrix(a, N, n_qubits) — Builds a $2^{n_work} \times 2^{n_work}$ permutation matrix for the map $z \mapsto (a \cdot z) \bmod P$ (identity for $z \geq P$). Used for both $g^{2^i}$ and $(y^{-1})^{2^j}$ controlled applications.
• _classical_post_processing(probs, g, y, P, n, N_size) — Multi-step classical routine:
  1. Sorts bitstrings by probability; skips entries below 0.02.
  2. Splits each bitstring into v_bin (first n bits, reg_a) and u_bin (last n bits, reg_b).
  3. Applies Fraction(u, N_size).limit_denominator(P) to get candidate (s_base, r_base).
  4. Searches multiples of r_base to find the true group order r where g^r ≡ 1 (mod P).
  5. Computes target = round(v * r / N_size) and solves $sx \equiv -\text{target} \pmod{r}$ via modular inverse, checking all solutions in the coset.
• _build_return_dict(is_success, circuit_path, filename, qc) — Packages the result dictionary returned by run(), including status ('ok'/'failed'), circuit_path, plot (list of file dicts), circuit, and the output fields Found x, Detected period r, and Computation time (s) merged from self.output.

Register address translation: The get_p(reg_slice) inline function inside run() translates named register slices into flat qubit indices by adding the appropriate offset (0 for reg_a, n_count for reg_b, 2*n_count for reg_work).

Data flow: (g, y, P) → register + oracle construction → qc.execute() → res_vec.calculate_state() → _classical_post_processing() → found_x → _build_return_dict().

Understanding the Key Quantum Components

Both reg_A and reg_B are placed in uniform superposition: $$\frac{1}{N}\sum_{a=0}^{N-1}\sum_{b=0}^{N-1}|a\rangle|b\rangle|1\rangle_{\text{work}}$$ where $N = 2^{n_{\text{count}}}$.

2. Modular Exponentiation (Oracle)

Phase 1: reg_A controls $\times g^{2^i}$ for each bit $i$: $$\sum_{a,b}|a\rangle|b\rangle|g^a \bmod P\rangle$$

Phase 2: reg_B controls $\times (y^{-1})^{2^j}$ for each bit $j$: $$\sum_{a,b}|a\rangle|b\rangle|g^a (y^{-1})^b \bmod P\rangle = \sum_{a,b}|a\rangle|b\rangle|g^{a-xb} \bmod P\rangle$$

3. Inverse QFT on Both Registers

The IQFT on both A and B converts the periodic structure into measurable peaks. The joint distribution peaks at $(u, v)$ where: $$\frac{u}{N} \approx \frac{s}{r}, \quad \frac{v}{N} \approx \frac{-sx}{r} \pmod{1}$$ for random $s \in {0, \ldots, r-1}$ (the random eigenstate index).

4. Phase Leakage into Momentum

The Dirichlet kernel describes how much of the measured probability concentrates near the true rational $s/r$: $$P(m) \approx \frac{\sin^2(\pi N \delta)}{N^2 \sin^2(\pi \delta)}, \quad \delta = m/N - s/r$$

5. Classical Post-Processing

From the measurement $(u, v)$:

1. Continued fractions on $u/N$ gives denominator $r$ (the group order) and numerator $s$.
2. Solve: $x \equiv -\text{round}(v \cdot r / N) \cdot s^{-1} \pmod{r}$ via modular inverse.

Theory-to-Code Mapping

Notes on encapsulation: The register address mapping (named registers → flat qubit indices) is handled by an inline get_p() closure inside run(), rather than a class-level method. This is because the three registers (reg_a, reg_b, reg_work) are passed as arguments to Circuit and the gate methods still require flat indices. The post-processing is entirely classical and self-contained in _classical_post_processing().

Mathematical Deep Dive

Group structure: The multiplicative group $(\mathbb{Z}/P\mathbb{Z})^*$ is cyclic of order $r = P-1$ when $P$ is prime and $g$ is a primitive root.

Eigenstate structure: The unitary $U_g|z\rangle = |gz \bmod P\rangle$ has eigenstates $|u_s\rangle = \frac{1}{\sqrt{r}}\sum_{k=0}^{r-1} e^{-2\pi i sk/r}|g^k\rangle$ with eigenvalues $e^{2\pi i s/r}$.

Two-dimensional QPE: Simultaneously doing QPE on $U_g$ (reg_A) and $U_y = U_g^x$ (reg_B) gives phase pair $(s/r, -sx/r)$, from which $x = -v \cdot (u)^{-1}$ can be recovered via modular arithmetic.

Complexity: $O(n^3)$ gates where $n = \lceil\log_2 P\rceil$. Classical best: sub-exponential $O(\exp(\sqrt{n\log n}))$ via index calculus.

Hands-On Example (UnitaryLab)

from unitarylab_algorithms import DiscreteLogAlgorithm

# Solve 3^x ≡ 6 (mod 7): answer is x=3 since 3^3=27=3*7+6
algo = DiscreteLogAlgorithm()
result = algo.run(g=3, y=6, P=7, backend='torch')
print(f"x = {result['Found x']}")                    # Expected: 3
print(f"Status: {result['status']}")                  # 'ok' on success, 'failed' on failure
print(f"Period r: {result['Detected period r']}")     # Detected group order
print(f"Time: {result['Computation time (s)']} s")    # Simulation time
print(f"Verify: 3^3 mod 7 = {pow(3, 3, 7)}")        # Should be 6
print(result.get('plot', []))  # [{"format": "txt", "filename": "..."}]

Implementing Your Own Version

Below is a skeleton that reconstructs the discrete-log quantum circuit at the component level, matching DiscreteLogAlgorithm.

# Simplified reconstruction — mirrors DiscreteLogAlgorithm.run()
import math
from fractions import Fraction
import numpy as np
from unitarylab import Circuit, Register
from unitarylab.library import IQFT

def modular_matrix(mult: int, P: int, n_work: int) -> np.ndarray:
    """Permutation matrix |x> -> |x * mult mod P>."""
    dim = 2**n_work
    U = np.zeros((dim, dim))
    for x in range(dim):
        y = (x * mult) % P if x < P else x
        U[y, x] = 1.0
    return U.astype(complex)

def build_dlp_circuit(g: int, y: int, P: int, n_count: int, n_work: int) -> Circuit:
    """
    Two-register QPE circuit for DLP g^x = y mod P.
    Register a (counts a): controlled g^{2^i} mod P
    Register b (counts b): controlled (y^{-1})^{2^j} mod P
    Work register: |1>, then |g^a * y^{-b}>
    Measurement: peak at (a, b) with a/r_g - b/r_y -> x via CRT
    """
    ra = Register('a', n_count)
    rb = Register('b', n_count)
    rw = Register('w', n_work)
    qc = Circuit(ra, rb, rw, name='DLP_circuit')

    # Offset helpers to get flat qubit indices
    a_off, b_off, w_off = 0, n_count, 2 * n_count

    # 1. Hadamard both counting registers; work register to |1>
    qc.h(range(a_off, a_off + n_count))
    qc.h(range(b_off, b_off + n_count))
    qc.x(w_off)   # |work> = |1>

    # 2. Controlled g^{2^i} on reg_a
    for i in range(n_count):
        mult = pow(g, 2**i, P)
        qc.unitary(modular_matrix(mult, P, n_work),
                   range(w_off, w_off + n_work), a_off + i, '1')

    # 3. Controlled (y^{-1})^{2^j} on reg_b
    y_inv = pow(y, -1, P)
    for j in range(n_count):
        mult = pow(y_inv, 2**j, P)
        qc.unitary(modular_matrix(mult, P, n_work),
                   range(w_off, w_off + n_work), b_off + j, '1')

    # 4. IQFT on both counting registers
    qc.append(IQFT(n_count), range(a_off, a_off + n_count))
    qc.append(IQFT(n_count), range(b_off, b_off + n_count))
    return qc

def solve_dlp(g: int, y: int, P: int, backend: str = 'torch', device: str = 'cpu'):
    """Quantum DLP: find x such that g^x = y mod P."""
    if math.gcd(g, P) != 1 or math.gcd(y, P) != 1:
        raise ValueError('g and y must be coprime to P')
    n_work  = P.bit_length()
    n_count = 2 * n_work

    qc = build_dlp_circuit(g, y, P, n_count, n_work)
    res_vec = qc.execute(backend=backend, device=device, dtype=np.complex128)
    probs_dict = res_vec.calculate_state(range(2 * n_count))

    # Sort by probability and post-process
    sorted_probs = sorted(probs_dict.items(), key=lambda item: item[1]['prob'], reverse=True)
    for bitstring, data in sorted_probs:
        if data['prob'] < 0.02:
            continue
        u_bin, v_bin = bitstring[n_count:], bitstring[:n_count]
        u, v = int(u_bin, 2), int(v_bin, 2)
        if u == 0:
            continue
        frac = Fraction(u, 2**n_count).limit_denominator(P)
        r = frac.denominator
        if r and pow(g, r, P) == 1:
            for x in range(r):
                if pow(g, x, P) == y:
                    return x
    return None

Component roles:

• modular_matrix — same pattern as Shor: $2^n \times 2^n$ permutation for $|x\rangle \to |x \cdot \text{mult} \bmod P\rangle$.
• build_dlp_circuit — two-register QPE: register $a$ accumulates powers of $g$, register $b$ accumulates powers of $y^{-1}$, both followed by IQFT.
• solve_dlp — measures both registers, extracts phases, and uses continued fractions to recover the periods that encode $x$.

Reference Implementation (Classiq)

Classiq provides a high-level QMOD implementation of the discrete logarithm algorithm. It uses @qfunc to define modular exponentiation, the discrete-log oracle, inverse qft, model creation, synthesis, and execution.

This reference is useful for understanding how the Shor-style discrete logarithm workflow can be expressed through automatic circuit synthesis.

Example A: Minimal Classiq Discrete Logarithm Model

from classiq import *
from classiq.qmod.symbolic import ceiling, log


@qfunc
def modular_exponentiation(N: CInt, a: CInt, x: QArray, pw: QArray):
    repeat(
        count=pw.len,
        iteration=lambda index: control(
            pw[index],
            lambda: modular_multiply_constant_inplace(
                N,
                a ** (2**index),
                x,
            ),
        ),
    )


@qfunc
def discrete_log_oracle(
    g_generator: CInt,
    x_element: CInt,
    N_modulus: CInt,
    alpha: QArray,
    beta: QArray,
    func_res: Output[QNum],
) -> None:
    allocate(ceiling(log(N_modulus, 2)), func_res)
    func_res ^= 1

    # Apply x^alpha mod N
    modular_exponentiation(N_modulus, x_element, func_res, alpha)

    # Apply g^beta mod N
    modular_exponentiation(N_modulus, g_generator, func_res, beta)


@qfunc
def discrete_log(
    g: CInt,
    x: CInt,
    N: CInt,
    order: CInt,
    alpha: Output[QArray],
    beta: Output[QArray],
    func_res: Output[QArray],
) -> None:
    reg_len = ceiling(log(order, 2))

    allocate(reg_len, alpha)
    allocate(reg_len, beta)

    hadamard_transform(alpha)
    hadamard_transform(beta)

    discrete_log_oracle(g, x, N, alpha, beta, func_res)

    invert(lambda: qft(alpha))
    invert(lambda: qft(beta))


MODULU_NUM = 5
G_GENERATOR = 3
X_LOG_ARG = 2
ORDER = MODULU_NUM - 1


@qfunc
def main(
    alpha: Output[QNum],
    beta: Output[QNum],
    func_res: Output[QNum],
) -> None:
    discrete_log(
        G_GENERATOR,
        X_LOG_ARG,
        MODULU_NUM,
        ORDER,
        alpha,
        beta,
        func_res,
    )


qmod = create_model(
    main,
    constraints=Constraints(max_width=13),
    preferences=Preferences(optimization_level=1),
    execution_preferences=ExecutionPreferences(num_shots=4000),
)

qprog = synthesize(qmod)
result = execute(qprog).result_value()
print(result)

Debugging Tips

1. Simulation is slow: Qubit count is $5\lfloor\log_2 P\rfloor$; exponential in bits. Keep $P$ small ($P=7$, $11$, $13$).
2. Found x is None: The continued fractions step failed to extract a valid period. Re-run (the algorithm is probabilistic) or increase $n_{\text{count}}$.
3. g and y not coprime to P: Will raise ValueError. Ensure $\gcd(g, P) = \gcd(y, P) = 1$.
4. Wrong answer: The quantum measurement gives the correct answer with high probability but not certainty. The classical post-processing verifies correctness; if wrong, retry.