Black-Box Optimization with Quantum Annealing and Ising Machines¶

This example program introduces a method of black-box optimization, which may be called Factorization Machines and Quantum Annealing (FMQA). You can find example programs of FMQA applications in various fields such as material search, fluid engineering, chemical plant and urban transportation.

Introduction to FMQA¶

Black-box optimization¶

FMQA is one of the black box optimization methods. Usually, in mathematical optimization, the objective is to find a set of decision variables $\boldsymbol{x}$ such that some objective function $y = f(\boldsymbol{x})$ is minimized (or maximized). Here, $\boldsymbol{x}$ is assumed to be a binary variable vector with size $d$ and each element taking the value 0 or 1.

$$ \begin{aligned} \mathrm{Minimize}&\,\,f(\boldsymbol{x}) \\ \mathrm{subject\,\,to\,\,}&\boldsymbol{x} \in [0,1]^d \end{aligned} $$

Here, if information about the objective function $f(\boldsymbol{x})$ (functional form, gradient, submodularity, convexity, etc.) is given, efficient optimization can be performed. For example, suppose the function $f(\boldsymbol{x})$ is known (and is quadratic in $\boldsymbol{x}$), as in some optimization problems shown in the Amplify demo tutorial. In such a case, $f(\boldsymbol{x})$ can be used as the objective function to perform the optimization directly as a quadratic unconstrained binary optimization (QUBO: Quadratic Unconstrained Binary Optimization) problem.

On the other hand, in the case of optimization to minimize (or maximize) values obtained by simulation or experiment for physical or social phenomena, the objective function $f(\boldsymbol{x})$ corresponds to simulation or experiment, and the function cannot be described explicitly. Mathematical optimization for such an unknown objective function $f(\boldsymbol{x})$ is called black-box optimization.

In addition, evaluating such an objective function (running simulations or experiments) is usually relatively expensive (in terms of time and money, etc). Therefore, even if the set of decision variables is finite, optimization by full search is generally difficult. Therefore, an optimization method with as few objective function evaluations as possible is required.

FMQA introduction¶

FMQA is a black-box optimization method that combines machine learning and quantum annealing. The method iterates through the cycles shown in the figure below, simultaneously searching for both a good approximation by a second-order polynomial of the black-box function and the input that yields the minimum value of that polynomial.

First, we compute a second-order polynomial that approximates the black box function using machine learning. Next, the optimization solver finds the input $x$ that minimizes the quadratic polynomial. The $x$ obtained by the Ising machine is then input to the black box function. Suppose the second-order polynomial obtained by machine learning approximates the black-box function well enough. In that case, we can expect that $\boldsymbol{x}$ outputs a small value when input to the black-box function. If not, we can still expect a better polynomial approximation of the black box function in the next training cycle by adding the data to the training data set and performing machine learning again.

The second-order polynomial approximation of the black-box function uses a machine-learning model called the Factorization Machine (FM), consisting of the following polynomial where $d$ is a constant representing the length of the input to the black-box function, $\boldsymbol{v}$, $\boldsymbol{w}$, and $w_0$ are the parameters of the model, and $k$ is a hyperparameter representing the size of the parameters.

$$ \begin{aligned} f(\boldsymbol{x} | \boldsymbol{w}, \boldsymbol{v}) &= w_0 + \langle \boldsymbol{w}, \boldsymbol{x}\rangle + \sum_{i=1}^d \sum_{j=i+1}^d \langle \boldsymbol{v}_i, \boldsymbol{v}_j \rangle x_i x_j \\ &=w_0 + \sum_{i=1}^d w_i x_i + \sum_{i=1}^d \sum_{j=i+1}^d \sum_{f=1}^k v_{if}v_{jf}x_ix_j \\ &=w_0 + \sum_{i=1}^d w_i x_i + \frac{1}{2}\sum_{f=1}^k\left(\left(\sum_{i=1}^d v_{i f} x_i\right)^2 - \sum_{i=1}^d v_{i f}^2 x_i^2\right) \end{aligned} $$

Using the FM as a machine learning model has the following advantages.

• Minimization by an optimization (QUBO) solver is possible because the model is a quadratic polynomial
• The computational complexity of the model's inference can be parameterized

The hyperparameter $k$ is a positive integer less than or equal to the length $d$ of the input to the black box function, which has the effect of adjusting the number of parameters in the FM model. When $k=d$, the model has the same degrees of freedom as the QUBO interaction term, while a smaller $k$ reduces the number of parameters and suppresses over-fitting.

FMQA procedure¶

In FMQA, we prepare the initial training data as follows, and then repeats the inference and optimization with the machine learning model.

• Preparation of initial training data ($N_0$ samples)

  1. Prepare $N_0$ input samples $\{\boldsymbol{x}_1, \boldsymbol{x}_2, \cdots, \boldsymbol{x}_{N_0}\}$ and the corresponding $N_0$ outputs $\{f(\boldsymbol{x}_1 ), f(\boldsymbol{x}_2), \cdots, \boldsymbol{x}_{N_0}\}$ as initial training data.

• FMQA optimization cycle ($N$ times)

  1. Train the FM model using the (most recent) training data and obtain the FM parameters $(\boldsymbol{v}, \boldsymbol{w})$.
  2. Estimate the input $\hat{\boldsymbol{x}}$ that minimizes the acquisition function $g(\boldsymbol{x})$ by using Amplify.
  3. Evaluate the objective function $f(\boldsymbol{x})$ with $\hat{\boldsymbol{x}}$ to obtain $\hat{y} = f(\hat{\boldsymbol{x}})$.
  4. Add $(\hat{\boldsymbol{x}}, \hat{y})$ to the training data.

  Repeat steps 1-4 above for $N$ times.

This example program shows how to run FMQA using PyTorch and the Amplify SDK. The part of the code that minimizes the learned model uses the Fixstars Amplify Annealing Engine (Amplify AE), a GPU-based annealing machine, instead of quantum annealing (QA).

The following is a example program that performs black-box optimization using FMQA.

Solver client configuration¶

First, we will configure the solver client of the solver we will use during the optimization cycles. In this example program, we will use Amplify AE.

Definition of black-box function¶

Now we define the black box function to be used for black box optimization. A black-box function is a function whose input is a 1-dimensional array of binary variables with values of 0 or 1 and whose output is a real number. However, input arrays consisting of integers and real numbers can also be considered. Examples of such use are presented here.

The black-box functions may simulate physical or social phenomena or return a value obtained from experiments for practical use. These are suitable for black-box optimization because they cannot be expressed in mathematical formulas, and their properties are not obvious.

However, in this tutorial, instead of simulation or experiment, we will prepare an simple algebraic function as a black box function. In the following, we define the make_blackbox_func function to create a "black-box" function that takes a NumPy 1-dimensional array with $d$ elements as input and real numbers as output.

The function created here is assumed to be a quadratic function, but from now on we will estimate and minimize the function, assuming that we do not know that the function created by make_blackbox_func is quadratic or any other properties.

Model training by machine learning¶

This section implements the part of FMQA that learns the optimal parameters (weights and bias) of the model by machine learning. It corresponds to the lower right part of the figure below, where the input is the training data and the output is the model.

!

First, the TorchFM class representing the model by the Factorization Machine is defined using PyTorch.

The following equation represents the Factorization Machine.

$$ \begin{aligned} f(\boldsymbol{x} | \boldsymbol{w}, \boldsymbol{v}) &= \underset{\color{red}{\mathtt{out\_linear}}}{\underline{ w_0 + \sum_{i=1}^d w_i x_i} } + \underset{\color{red}{\mathtt{out\_quadratic}}}{\underline{\frac{1}{2} \left[\underset{\color{red}{\mathtt{out\_1}}}{\underline{ \sum_{f=1}^k\left(\sum_{i=1}^d v_{i f} x_i\right)^2 }} - \underset{\color{red}{\mathtt{out\_2}}}{\underline{ \sum_{f=1}^k\sum_{i=1}^d v_{i f}^2 x_i^2 }} \right] }} \end{aligned} $$

The input $x$ of this model is a vector of the length $d$ as the input to the black-box function, with the following three parameters.

• $v$: 2-dimensional array of $d\times k$.
• $w$: 1D vector of length $d$.
• $w_0$: scalar

The only hyperparameter is $k$, which is given as a positive integer less than or equal to $d$.

The TorchFM class defined below inherits from torch.nn.Module and is constructed from an input vector $x$ of size $d$ and a hyperparameter $k$. The hyperparameter $k$ controls the number of parameters in the model; the larger the hyperparameter, the more parameters are, the more accurate the model becomes, but also the more prone to over-fitting.

The TorchFM class has attributes $v$, $w$, and $w_0$ of the model parameters and updates these parameters as the training proceeds. According to the above formula, the forward method also outputs an estimate of $y$ from the input $x$. Since the parameters $v$, $w$, and $w_0$ are needed to construct a QUBO model later, we also define a function get_parameters to output these parameters.

Next, we define a function train() to perform machine learning of the model based on the TorchFM class. The input is the training data $x, y$ and an instance of the TorchFM model. By calling the train() function, the TorchFM parameters are trained and optimized.

In general machine learning, we split the data into training and validation data, optimize the parameters using the training data, and validate the model during training using the validation data. The model is validated at each epoch, and the parameters at the epoch with the highest prediction accuracy for the validation data are saved and used as the model after learning.

Solving the model with Amplify¶

Next, we implement the anneal function to minimize the inferred machine learning model, corresponding to the lower left part of the FMQA cycle diagram (reproduced below), where the input is the model TorchFM class after training, and the output is a vector $x$ that minimizes the model.

Solve the following optimization problem for the class TorchFM model learned earlier to find the input $x$ such that the inferred model is minimized.

$$ \underset{x}{\mathrm{argmin}} \quad \underset{\color{red}{\mathtt{out\_linear}}}{\underline{ w_0 + \sum_{i=1}^d w_i x_i} } + \underset{\color{red}{\mathtt{out\_quadratic}}}{\underline{\frac{1}{2} \left[\underset{\color{red}{\mathtt{out\_1}}}{\underline{ \sum_{f=1}^k\left(\sum_{i=1}^d v_{i f} x_i\right)^2 }} - \underset{\color{red}{\mathtt{out\_2}}}{\underline{ \sum_{f=1}^k\sum_{i=1}^d v_{i f}^2 x_i^2 }} \right] }} $$

The decision variable in this optimization problem is $x$. It is a one-dimensional binary variable vector of length $d$, just like the input vector to the black-box function. Also, $v$, $w$, and $w_0$, parameters in the learning phase, are constants here.

We define the anneal function that performs the optimization using Amplify AE for the given model as follows. The anneal function uses the VariableGenerator to create a 1D binary vector x of length $d$, and then uses the binary variable array $x$ and $v$, $w$, and $w_0$ obtained from the TorchFM class to create a Factorization Machine The objective function to be optimized is created according to the formula of Factorization Machine.

After building the optimization model, we will perform minimization of the objective function using the created QUBO model and the solver client FixstarsClient defined earlier.

Performing FMQA¶

Now that we have defined the train function for machine learning and the anneal function for optimization, which are the core of FMQA, we can execute FMQA using these functions.

First, we create a black-box function for black-box optimization as follows. This function takes a NumPy one-dimensional vector of length $d = 100$ consisting of $0$ or $1$ and returns a float.

Creating initial training data¶

Next, the black-box function $y = f(\boldsymbol{x})$ is evaluated on the input vector $\boldsymbol{x}$ to create $N_0$ initial training data. The black-box function corresponds to the results of a simulation or experiment, so you can create these data based on previous experiments or simulations.

Since we have prepared an appropriate function as a black-box function for simulation, we define the init_training_data function to create the initial training data using a random $N_0$ input vector $x$ as follows. This function takes the black-box function and the number of initial training data $N_0$ and returns $N_0$ input vectors $\boldsymbol{x}$ and corresponding output $y$ as initial training data.

Now, we have $N_0$-pair of initial training data.

Running an FMQA cycle¶

Using the $N_0$ pairs of initial training data created above, we will run the FMQA cycle according to the following diagram.

We perform the following operations per FMQA cycle.

1. Train the model
   • Construct a model of class TorchFM and train it by calling the train function on the model with initial training data x = x_init, y = y_init.
2. Minimize the model
   • Minimize the model by passing a trained model of class TorchFM to the anneal function to get $\hat{x}$ that minimizes the model.
   • If $\hat{x}$ is already included in the training data x, change part of $\hat{x}$ to avoid duplicate training data.
3. Evaluation of black-box function
   • Input $\hat{x}$ to the black-box function and get output $\hat{y}$.
   • Add $\hat{x}$ and $\hat{y}$ to the traiing data x and y, respectively.

An example implementation of the above is as follows: Since the execution takes a few minutes of computation time, an example output is shown in Execution example.

After running the above cell, x and y have input-output pairs of black-box function evaluations. The number of pairs is $N_0 + N = 70$.

The following gives the input that gives the minimum value of the black box function evaluation and its value.

Plotting the evolution of the evaluated values¶

Below are $N_0$ initial training data and the evolution of the optimized evaluation values of the black-box function over $N$ FMQA cycles. The initial training data is shown in blue, and the evaluation values obtained by the FMQA cycles are shown in red.

Execution example¶

The output of a typical plot obtained by running this sample code is shown below. The optimization uses a QUBO solver, so the results will vary from run to run. However, you can see that the evaluated value of the black-box function decreases with the number of FMQA cycles.

4. References¶

The present black-box optimization method that combines quantum annealing and Ising machines with machine learning is called FMQA, which has been originally proposed as FMQA in the following research.

• K. Kitai, J. Guo, S. Ju, S. Tanaka, K. Tsuda, J. Shiomi, and R. Tamura, "Designing metamaterials with quantum annealing and factorization machines", Physical Review Research 2, 013319 (2020).

In this study, the search for "metamaterials" is carried out using FMQA, which also have shown superior performance compared to Bayesian optimization, a conventional black-box optimization method.

In the following study, the same black-box optimization method is also applied to the design of photonic crystals.

• T. Inoue, Y. Seki, S. Tanaka, N. Togawa, K. Ishizaki, and S. Noda, "Towards optimization of photonic-crystal surface-emitting lasers via quantum annealing," Opt. Express 30, 43503-43512 (2022).

These studies suggest that this optimization method (FMQA), based on FM and combinatorial optimization, may have general applicability in black-box optimization problems in various fields. In Fixstars Amplify, there are several examples of such black-box optimization in the areas of chemical reaction, fluid dynamics, as well as material search, as follows:

• Black-Box Optimization Exploration of Model Superconducting Materials
• Black-Box Optimization of Operating Condition in a Chemical Reactor
• Black-Box Optimization of Airfoil Geometry by Fluid Flow Simulation