Prime Factorization

A central problem in Number Theory and cryptography is to factor $R \in \mathbb{N}$, which is known to be the product of two distinct prime numbers $p, q \in \mathbb{N}$. Shor's Algorithm, designed to address such task is often regarded as one of the major theoretical landmarks in Quantum Computing, being responsible for driving increasingly greater interest to the area.

QUBO Formulation

A naïve approach to model this problem can be stated as a quadratically-constrained integer program:

\[\begin{array}{rl} \text{s.t.} & p \times q = R \\ & p, q \ge 0 \\ & p, q \in \mathbb{Z} \end{array}\]

This is a constraint satisfaction problem (no objective function) where we seek integer values $p$ and $q$ that multiply to give $R$.

Variable Bounds

From the definition and the basics of number theory, we are able to retrieve a few assumptions about the problem's variables:

• $p$ and $q$ are bounded to the interval $\left[1, R\right]$
• Moreover, it is fine to assume $1 < p \le \sqrt{R} \le q \le R \div 2$.

These bounds are crucial for the QUBO encoding, as they determine how many binary variables are needed to represent $p$ and $q$.

Implementation

using JuMP
using ToQUBO
using PySA

function factor(R::Integer; optimizer = PySA.Optimizer)
    return factor(identity, R; optimizer)
end

function factor(config!::Function, R::Integer; optimizer = PySA.Optimizer)
    model = Model(() -> ToQUBO.Optimizer(optimizer))

    @variable(model,  1 <= p <= √R, Int)
    @variable(model, √R <= q <= R ÷ 2, Int)

    @constraint(model, p * q == R)

    config!(model)

    optimize!(model)

    p = round(Int, value(model[:p]))
    q = round(Int, value(model[:q]))

    return (p, q)
end

factor (generic function with 2 methods)

Solving

p, q = factor(15) do model
    set_optimizer_attribute(model, "n_reads", 100)
    set_optimizer_attribute(model, "n_sweeps", 200)
end

print("$p, $q")

  0%|          | 0/100 [00:00<?, ?it/s]/home/runner/work/ToQUBO.jl/ToQUBO.jl/docs/.CondaPkg/.pixi/envs/default/lib/python3.11/site-packages/pysa/simulation.py:157: NumbaPerformanceWarning: np.dot() is faster on contiguous arrays, called on (Array(float64, 1, 'A', False, aligned=True), Array(float64, 1, 'C', False, aligned=True))
  energies[k] += sweep(update_spin, couplings, local_fields,

  1%|          | 1/100 [00:07<12:58,  7.86s/it]
100%|██████████| 100/100 [00:07<00:00, 12.57it/s]
3, 5

Scaling Considerations

While this approach works well for small numbers, the number of binary variables required grows with the size of the factors. For a number $R$ with $n$-bit factors, approximately $O(n)$ binary variables are needed for each factor. The quadratic constraint $p \times q = R$ then introduces $O(n^2)$ terms in the QUBO, which determines the problem's complexity on quantum hardware.