Multiplier Simulation and Factorization
Multiplication of two integers can be performed using additions. In this section, we design a multiplier for two 4-bit integers using full adders. The figure below shows how two $x_3x_2x_1x_0$ and $y_3y_2y_1y_0$ are multiplied to obtain an 8-bit integer $z_7z_6z_5z_4z_3z_2z_1z_0$. In this figure, $p_{i,j}=x_iy_j$ ($0\leq i,j\leq 3$) and these partial products are summed to compute the final 8-bit result.
We use a 4-bit ripple-carry adder that computes the sum of two 4-bit integers $a_3a_2a_1a_0$ and $b_3b_2b_1b_0$ producing the 5-bit sum $z_4z_3z_2z_1z_0$. It consists of four full adders connected by a 5-bit carry wire $c_4c_3c_2c_1c_0$ that propagates carries.
A 4-bit multiplier can be constructed using three 4-bit adders. They are connected by wires $c_{i,j}$ ($0\leq i\leq 2, 0\leq j\leq 3$) to propagate intermediate sum bits, as shown below:
QUBO formulation for multiplier
We will show QUBO formulation for simulating the N-bit multiplier. To do this, we implement functions that construct a full adder, an adder, and a multiplier.
Full adder
The following QUBO expression simulates a full adder with three input bits a, b, and i, and two output bits: carry-out o and sum s:
qbpp::Expr fa(const qbpp::Expr& a, const qbpp::Expr& b, const qbpp::Expr& i,
const qbpp::Expr& o, const qbpp::Expr& s) {
return (a + b + i) - (2 * o + s) == 0;
}
The function fa returns an expression that enforces consistency between the input and output bits of a full adder.
Adder
Assume that arrays a, b, and s of qbpp::Expr objects represent integers. We assume that a and b each have N elements representing N-bit integers, while s has N + 1 elements representing an (N + 1)-bit integer. The following function adder returns a QUBO expression whose minimum value is 0 if and only if a + b == s holds:
qbpp::Expr adder(const qbpp::ArrayBase& a,
const qbpp::ArrayBase& b,
const qbpp::ArrayBase& s) {
auto N = a.size();
auto c = qbpp::var("_c", N + 1);
auto f = qbpp::toExpr(0);
for (size_t j = 0; j < N; ++j) {
f += fa(qbpp::Expr(a[j]), qbpp::Expr(b[j]), qbpp::Expr(c[j]), qbpp::Expr(c[j + 1]), qbpp::Expr(s[j]));
}
f.replace({{qbpp::Var(c[0]), 0}, {qbpp::Var(c[N]), qbpp::Expr(s[N])}});
return f;
}
In this function, c is an array of N + 1 variables used to connect the carry-out and carry-in signals of the fa blocks, forming an N-bit ripple-carry adder.
Multiplier
Assume that arrays x, y, and z of qbpp::Expr represent integers. We assume that x and y each have N elements and that z has 2 * N elements. The following function multiplier returns a QUBO expression whose minimum value is 0 if and only if x * y == z holds.
qbpp::Expr multiplier(const qbpp::ArrayBase& x,
const qbpp::ArrayBase& y,
const qbpp::ArrayBase& z) {
auto N = x.size();
auto c = qbpp::var("c", N - 1, N + 1);
auto f = qbpp::toExpr(0);
for (size_t i = 0; i < N - 1; ++i) {
auto b = qbpp::expr(N);
for (size_t j = 0; j < N; ++j) {
b.at(j) = qbpp::Expr(x[i + 1]) * qbpp::Expr(y[j]);
}
auto a = qbpp::expr(N);
if (i == 0) {
for (size_t j = 0; j < N - 1; ++j) {
a.at(j) = qbpp::Expr(x[0]) * qbpp::Expr(y[j + 1]);
}
a.at(N - 1) = 0;
} else {
for (size_t j = 0; j < N; ++j) {
a.at(j) = qbpp::Expr(c[i - 1][j + 1]);
}
}
auto s = qbpp::expr(N + 1);
for (size_t j = 0; j < N + 1; ++j) {
s.at(j) = qbpp::Expr(c[i][j]);
}
f += adder(a, b, s);
}
f += qbpp::Expr(z[0]) - qbpp::Expr(x[0]) * qbpp::Expr(y[0]) == 0;
qbpp::MapList ml;
for (size_t i = 0; i < N - 2; ++i) {
ml.push_back({qbpp::Var(c[i][0]), qbpp::Expr(z[i + 1])});
}
for (size_t i = 0; i < N + 1; ++i) {
ml.push_back({qbpp::Var(c[N - 2][i]), qbpp::Expr(z[N + i - 1])});
}
return f.replace(ml).simplify_as_binary();
}
This function uses an (N−1)×(N+1) matrix c of qbpp::Var objects to connect the N−1 adders of N bits. Since each bit of z corresponds to one element of c, their correspondence is defined in ml, and the replacements are performed using replace().
QUBO++ program for factorization
Using the function multiplier, we can factor a composite integer into two factors. The following program constructs a 4-bit multiplier with
x: 4 binary variables,y: 4 binary variables,z: an array of constants{1, 1, 1, 1, 0, 0, 0, 1}, representing the 8-bit integer10001111(143), and stores the resulting expression inf:
#include <qbpp/qbpp.hpp>
#include <qbpp/easy_solver.hpp>
qbpp::Expr fa(const qbpp::Expr& a, const qbpp::Expr& b, const qbpp::Expr& i,
const qbpp::Expr& o, const qbpp::Expr& s) {
return (a + b + i) - (2 * o + s) == 0;
}
qbpp::Expr adder(const qbpp::ArrayBase& a,
const qbpp::ArrayBase& b,
const qbpp::ArrayBase& s) {
auto N = a.size();
auto c = qbpp::var("_c", N + 1);
auto f = qbpp::toExpr(0);
for (size_t j = 0; j < N; ++j) {
f += fa(qbpp::Expr(a[j]), qbpp::Expr(b[j]), qbpp::Expr(c[j]), qbpp::Expr(c[j + 1]), qbpp::Expr(s[j]));
}
return f.replace({{qbpp::Var(c[0]), 0}, {qbpp::Var(c[N]), qbpp::Expr(s[N])}});
}
qbpp::Expr multiplier(const qbpp::ArrayBase& x,
const qbpp::ArrayBase& y,
const qbpp::ArrayBase& z) {
auto N = x.size();
auto c = qbpp::var("c", N - 1, N + 1);
auto f = qbpp::toExpr(0);
for (size_t i = 0; i < N - 1; ++i) {
auto b = qbpp::expr(N);
for (size_t j = 0; j < N; ++j) {
b.at(j) = qbpp::Expr(x[i + 1]) * qbpp::Expr(y[j]);
}
auto a = qbpp::expr(N);
if (i == 0) {
for (size_t j = 0; j < N - 1; ++j) {
a.at(j) = qbpp::Expr(x[0]) * qbpp::Expr(y[j + 1]);
}
a.at(N - 1) = 0;
} else {
for (size_t j = 0; j < N; ++j) {
a.at(j) = qbpp::Expr(c[i - 1][j + 1]);
}
}
auto s = qbpp::expr(N + 1);
for (size_t j = 0; j < N + 1; ++j) {
s.at(j) = qbpp::Expr(c[i][j]);
}
f += adder(a, b, s);
}
f += qbpp::Expr(z[0]) - qbpp::Expr(x[0]) * qbpp::Expr(y[0]) == 0;
qbpp::MapList ml;
for (size_t i = 0; i < N - 2; ++i) {
ml.push_back({qbpp::Var(c[i][0]), qbpp::Expr(z[i + 1])});
}
for (size_t i = 0; i < N + 1; ++i) {
ml.push_back({qbpp::Var(c[N - 2][i]), qbpp::Expr(z[N + i - 1])});
}
return f.replace(ml).simplify_as_binary();
}
int main() {
auto x = qbpp::var("x", 4);
auto y = qbpp::var("y", 4);
auto z = qbpp::int_array({1, 1, 1, 1, 0, 0, 0, 1});
auto f = multiplier(x, y, z).simplify_as_binary();
auto solver = qbpp::easy_solver::EasySolver(f);
auto sol = solver.search({{"target_energy", 0}});
for (size_t i = x.size(); i > 0; --i) {
std::cout << sol(x[i - 1]);
}
std::cout << " * ";
for (size_t i = y.size(); i > 0; --i) {
std::cout << sol(y[i - 1]);
}
std::cout << " = ";
for (size_t i = z.size(); i > 0; --i) {
std::cout << z[i - 1];
}
std::cout << std::endl;
}
The Easy Solver is executed on f, and the obtained solution is stored in sol. The resulting values of x and y are printed as:
1011 * 1101 = 10001111
This output indicates $11\times 13 = 143$, demonstrating the factorization result.