qbpp_multiplier.hpp

Go to the documentation of this file.

 
#ifndef QBPP_MULTIPLIER_HPP
#define QBPP_MULTIPLIER_HPP
 
#include "qbpp.hpp"
 
namespace qbpp {
namespace factorization {
 
inline Expr and_gate(const Expr &a, const Expr &b, Var x, MapList &opt_sol) {
  opt_sol.push_back({x, eval(a, opt_sol) & eval(b, opt_sol)});
  return a * b - 2 * a * x - 2 * b * x + 3 * x;
}
 
inline Expr half_adder(const Expr &x, const Expr &y, Var s, Var c,
                       MapList &opt_sol) {
  auto val_x = eval(x, opt_sol);
  auto val_y = eval(y, opt_sol);
  opt_sol.push_back({s, val_x ^ val_y});
  opt_sol.push_back({c, val_x & val_y});
  return s + x + y + 4 * c - 4 * c * x - 4 * c * y - 2 * s * x - 2 * s * y +
         2 * x * y + 4 * c * s;
}
 
inline Expr full_adder(const Expr &x, const Expr &y, const Expr &z, Var s,
                       Var c, MapList &opt_sol) {
  auto val_x = eval(x, opt_sol);
  auto val_y = eval(y, opt_sol);
  auto val_z = eval(z, opt_sol);
  opt_sol.push_back({s, val_x ^ val_y ^ val_z});
  opt_sol.push_back({c, (val_x & val_y) | (val_y & val_z) | (val_z & val_x)});
  return s + x + y + z + 4 * c - 4 * c * x - 4 * c * y - 4 * c * z - 2 * s * x -
         2 * s * y - 2 * s * z + 2 * x * y + 2 * x * z + 2 * y * z + 4 * c * s;
}
 
template <typename T, typename U>
Expr adder(const qbpp::Vector<T> &x, const qbpp::Vector<U> &y,
           const qbpp::Vector<Var> &s, const qbpp::Vector<Var> &c,
           MapList &opt_sol) {
  size_t n = x.size();
  if (n != y.size() || n != s.size() || n != c.size()) {
    throw std::runtime_error("Invalid input size");
  }
  Expr expr = half_adder(x[0], y[0], s[0], c[0], opt_sol);
  for (size_t i = 1; i < x.size(); i++) {
    expr += full_adder(x[i], y[i], c[i - 1], s[i], c[i], opt_sol);
  }
  return expr;
}
 
template <typename T, typename U>
Expr multiplier(const qbpp::Vector<T> &x, const qbpp::Vector<U> &y,
                const qbpp::Vector<Var> &z, MapList &opt_sol) {
  int Xsize = static_cast<int>(x.size());
  int Ysize = static_cast<int>(y.size());
  int Zsize = static_cast<int>(z.size());
  if (Xsize + Ysize != Zsize) {
    throw std::runtime_error("Invalid input size");
  }
 
  // The resulting QUBO expression is stored in expr.
  Expr expr;
 
  // Creates variable array p of size Ysize * Xsize
  auto p = var("p", Ysize, Xsize);
  // Takes minimum value 0 if p[i][j] = y[i] * x[j].
  for (int i = 0; i < Ysize; i++)
    for (int j = 0; j < Xsize; j++) {
      expr += and_gate(y[i], x[j], p[i][j], opt_sol);
    }
  // Creates variable array s of size Ysize - 1 * Xsize
  // s is the output s of a half/full adder in an Xsize-bit adder
  // Each s[i] is a vector of the output s of the i-th adder.
  auto s = var("s", Ysize - 1, Xsize);
  // Creates variable array c of size Ysize * Xsize
  // c is the output c of a half/full adder in an Xsize-bit adder
  // Each c[i] is a vector of the output c of the i-th adder.
  auto c = var("c", Ysize, Xsize);
 
  // Temporary vector P of variables to be given to the 0-th adder
  qbpp::Vector<Expr> P(p[0].begin() + 1, p[0].end());
  P.emplace_back(0);
 
  // QUBO expression for the 0-th adder
  expr += adder(P, p[1], s[0], c[0], opt_sol);
  // QUBO expression for the i-th adder (1 <= i < Ysize - 1)
  for (int i = 1; i < Ysize - 1; ++i) {
    qbpp::Vector<Expr> S(s[i - 1].begin() + 1, s[i - 1].end());
    S.emplace_back(c[i - 1][Xsize - 1]);
    expr += adder(p[i + 1], S, s[i], c[i], opt_sol);
  }
 
  // p[0][0] is replaced by z[0]
  MapList map = {{p[0][0], z[0]}};
  // s[1][0], s[2][0],..., s[Ysize - 2][0] are replaced by z[1],
  // z[2],...,z[Ysize - 2]
  for (int i = 1; i < Ysize; ++i) {
    map.push_back({s[i - 1][0], z[i]});
  }
  // s[Ysize - 2][1], s[Ysize - 2][2],..., s[Ysize - 2][Xsize - 1] are replaced
  // by z[Ysize - 1],  z[Ysize],...,z[Xsize + Ysize - 2]
  for (int j = 1; j < Xsize; ++j) {
    map.push_back({s[Ysize - 2][j], z[j + Ysize - 1]});
  }
  // c[Ysize - 2][Xsize - 1] is replaced by z[Xsize + Ysize - 1]
  map.push_back({c[Ysize - 2][Xsize - 1], z[Xsize + Ysize - 1]});
 
  // Replacement is performed and then simplified as a binary expression.
  expr.replace(map).simplify_as_binary();
 
  return expr;
}
}  // namespace factorization
}  // namespace qbpp
 
#endif