Optimization method and information processing device
By relaxing binary variables with k-hot constraints to continuous values and performing parallel updates, the method efficiently addresses combinatorial optimization problems with k-hot constraints, enhancing the speed of Markov Chain Monte Carlo and simulated annealing processes.
Patent Information
- Application Number
- JP2023002594
- Authority / Receiving Office
- JP · JP
- Patent Type
- Patents
- Current Assignee / Owner
- Filing Date
- 2023-01-11
- Publication Date
- 2026-03-05
- Estimated Expiration
- 2043-01-11
AI Technical Summary
Existing methods struggle to efficiently solve combinatorial optimization problems with k-hot constraints on binary variables, such as those encountered in work shift optimization and delivery planning, due to the interdependence of variables and constraints, which hinders the speed of Markov Chain Monte Carlo and simulated annealing processes.
The method involves relaxing binary variables with k-hot constraints to continuous values and adding a correction value to maintain a constant sum, allowing for state transitions that satisfy the constraints, and utilizing an information processing device with an energy calculation unit and connection strength calculation unit to perform parallel updates.
This approach enables efficient solving of quadratic programming problems with k-hot constraints by allowing independent updates of variables, thereby speeding up Markov Chain Monte Carlo and simulated annealing processes.
Smart Images

Figure 0007824897000016 
Figure 0007824897000017 
Figure 0007824897000018
Abstract
Description
[Technical Field]
[0001] The present invention relates to an optimization method and an information processing device. [Background technology]
[0002] Patent Document 1 describes a method for efficiently searching for ground states based on simulated annealing by parallelizing energy calculations when k-hot constraints are imposed on variables of an Ising model with a quadratic energy function.
[0003] Patent Document 2 describes a method for efficiently searching for ground states based on simulated annealing by converting the interaction relationship of an Ising model with a quadratic energy function into a complete bipartite graph structure.
[0004] Patent Document 3 describes a method for efficiently searching for a ground state based on simulated annealing by converting the interaction relationship into a complete bipartite graph structure for a mixed binary quadratic programming problem that includes the ground state search problem of the Ising model.
[0005] Non-Patent Document 1 describes a method of random sampling from a truncated normal distribution whose domain is bounded on both the upper and lower sides. [Prior art documents] [Patent documents]
[0006] [Patent Document 1] Japanese Patent Publication No. 2020-64535 [Patent Document 2] WO2019 / 216277A1 [Patent Document 3] WO2021 / 220445A1 [Non-patent literature]
[0007] [Non-Patent Document 1] ZI Botev, “The normal law under linear restrictions: simulation and estimation via minimax tilting”, Journal of the Royal Statistical Society: Series B (2016). Summary of the Invention [Problem to be solved by the invention]
[0008] Many physical and social phenomena can be expressed using interaction models. An interaction model is defined by multiple nodes that make up the model, nonlinear coefficients (or interaction coefficients) between the nodes, and, if necessary, linear coefficients (or bias coefficients) that act on each node. In the fields of physics and social science, various models, including the Ising model, have been proposed, and all of them can be interpreted as a form of interaction model.
[0009] Combinatorial optimization problems, which involve finding a combination of variables that minimizes (or maximizes) an evaluation function, can also be expressed as a problem of searching for the state (ground state) of the node that minimizes the evaluation function by assigning variables to each node. Furthermore, real-world combinatorial optimization problems often involve conditions that the solution must satisfy, and these conditions become constraints imposed between nodes. Therefore, in problem solving, it is important to efficiently find a state that minimizes the evaluation function while satisfying the constraints.
[0010] Among these constraints, one that frequently appears in real-world combinatorial optimization problems is the k-hot constraint. A k-hot constraint is a constraint that requires that k items (k is a positive integer) in a set of binary variables with two states, on / off, must be on. For example, a k-hot constraint appears when assigning a specified number of items or events to a person in charge, such as in work shift optimization problems or delivery planning problems.
[0011] The present invention has been made in view of the above background, and aims to efficiently solve a quadratic programming problem having k-hot constraints (k is a positive integer) on binary variables. [Means for solving the problem]
[0012] A preferred aspect of the present invention is an optimization method for solving a quadratic programming problem in which one or more independent k-hot constraints (k is a positive natural number) are imposed on binary variables using an information processing device having a processor, a memory device, an input device, and an output device, in which the information processing device adds a correction value to the nonlinear coefficient of the binary variables to which the k-hot constraints are imposed to relax the binary variables to continuous values, and the information processing device performs state transitions so that the sum of the set of continuous-valued variables to which the k-hot constraints are imposed becomes constant, thereby performing a solution search while satisfying the k-hot constraints.
[0013] Another preferred aspect of the present invention is an information processing device that includes a processor, a memory device, an input device, an output device, and an arithmetic unit, and that solves a quadratic programming problem in which one or more independent k-hot constraints (k is a positive natural number) are imposed on binary variables.The information processing device includes an energy calculation execution unit that performs continuous value relaxation on the binary variables to which the k-hot constraints are imposed, and a connection strength calculation unit that calculates a correction value to be added to the nonlinear coefficient of the binary variables to which the k-hot constraints are imposed.Under the control of the energy calculation execution unit, the arithmetic unit performs state transitions so that the sum of the set of continuous value variables to which the k-hot constraints are imposed is constant, thereby performing a solution search while satisfying the k-hot constraints.
[0014] Other problems and solutions disclosed in the present application will be made clear in the detailed description and drawings. [Effects of the Invention]
[0015] According to the present invention, it is possible to efficiently solve a quadratic programming problem with k-hot constraints (k is a positive integer) on binary variables. Other problems, configurations, and advantages will become clear from the following description of the embodiments of the invention. [Brief explanation of the drawings]
[0016] [Figure 1] Schematic diagram of the energy landscape of the interaction model. [Figure 2] This is a conceptual diagram in which the interaction relationships between variables in the interaction model are represented as a complete graph, and binary variables are relaxed to continuous values. [Figure 3] This is a conceptual diagram of converting an interaction model represented by a complete graph into an alternative model in which the adjacency matrix between nodes has a complete bipartite graph structure in order to alternatively perform a ground state search for the interaction model represented by a complete graph. [Figure 4] FIG. 4 is a diagram illustrating a method for updating variables in the interaction model having the complete bipartite graph structure of FIG. 3. [Figure 5A] FIG. 10 is a graph showing the range of motion of xi when the sum of the pair of variables is kept constant. [Figure 5B] FIG. 5B is a graphical representation of the probability density function followed by the variables over the range of motion of FIG. 5A. [Figure 6] FIG. 1 is a block diagram showing a schematic configuration of an information processing device. [Figure 7] FIG. 2 is a block diagram of an arithmetic circuit. [Figure 8] FIG. 2 is a functional block diagram showing main functions of the information processing device. [Figure 9] 10 is a flowchart illustrating a ground state search process. [Figure 10] FIG. 2 is a block diagram illustrating a detailed configuration example of a calculation device. DETAILED DESCRIPTION OF THE INVENTION
[0017] Hereinafter, embodiments will be described in detail with reference to the drawings. In the following description, identical or similar components will be designated by common reference numerals, and duplicate descriptions will be omitted. Furthermore, when there are multiple elements having the same or similar functions, different subscripts will be used for the same reference numerals. Furthermore, when there is no need to distinguish between multiple elements, the subscripts will be omitted.
[0018] The designations "first," "second," "third," etc. in this specification are used to identify components and do not necessarily limit the number, order, or content thereof. Furthermore, numbers used to identify components are used in different contexts, and numbers used in one context do not necessarily indicate the same configuration in another context. Furthermore, this does not prevent a component identified by a certain number from also serving the function of a component identified by another number.
[0019] In order to facilitate understanding of the invention, the position, size, shape, range, etc. of each component shown in the drawings etc. may not represent the actual position, size, shape, range, etc. Therefore, the present invention is not necessarily limited to the position, size, shape, range, etc. disclosed in the drawings etc.
[0020] All publications, patents, and patent applications cited herein are incorporated by reference in their entirety.
[0021] As used herein, elements referred to in the singular are intended to include the plural unless the context clearly indicates otherwise.
[0022] In this embodiment, an optimization method for solving a mixed binary quadratic programming problem having one or more independent k-hot constraints on binary variables using an information processing device is described. In this method, all binary variables are relaxed to continuous values after correction is made to nonlinear coefficients for variables with k-hot constraints. The ground state search for this continuous-value relaxed interaction model with k-hot constraints is replaced with a ground state search for an interaction model whose adjacency matrix has a bipartite graph structure and in which k-hot constraints are imposed between variables. The information processing device includes an execution unit that searches for the ground state of the interaction model based on an energy function of the interaction model and information on the variables.
[0023] Other problems and solutions disclosed in the present application will be made clear in the detailed description and drawings.
[0024] First, we will explain the interaction model. The interaction model is defined by multiple nodes that make up the model, nonlinear coefficients acting between the nodes, and linear coefficients acting on each node. Here, we define a variable s corresponding to each node i (i = 1 to N). i is a continuous variable s- i ∈[-1,1] or binary variable s- i ∈{-1,1}. Based on the nonlinear and linear coefficients of the interaction model, the energy function H(s) (also called Hamiltonian) is defined. In particular, when the interaction relationship is defined between two nodes, the energy function is the following quadratic equation:
[0025]
number
[0026] It can also be thought of as an energy function consisting of the first term in Equation 1 representing the interaction between nodes, and the second term representing the bias for the node. Generally, an interaction model is expressed as an undirected graph, and the nonlinear terms are specified by the set of variables contained therein. Therefore, the matrix J in Equation 1 is a real symmetric matrix.
[0027] The above interaction model is a concept that includes the Ising model. In the interaction model, all variables are treated as binary variables s i When we restrict it to ∈{-1,1}, we get the Ising model, which is used, for example, as a lattice model to describe magnetic materials in statistical mechanics, with +1 / -1 corresponding to up / down spin.
[0028] The ground state search for the interaction model is a combinatorial optimization problem for finding a variable array s that minimizes the above energy function. In this embodiment, the ground state search for the interaction model is performed using Markov Chain Monte Carlo methods (hereinafter referred to as MCMC).
[0029] Figure 1 is a conceptual diagram of the energy landscape of an interaction model. The horizontal axis of the graph is the variable array, and the vertical axis is the energy of the system. In a state transition according to MCMC, a stochastic transition from the current state s to another state is repeated. Methods for sequentially executing stochastic state transitions include Gibbs sampling (or the heat bath method) and the Metropolis method. Let p(s) be the probability distribution function that state s follows, and s be the state generated at the i-th step. (i) ={s1 (i) , …, s N (i) In particular, in Gibbs sampling, the state of the i+1th iteration is given by the conditional probability
[0030]
number
[0031] As an example of such a sequential transition, in the case of Figure 1, -1 in state A is inverted to 1 to become state B, and then one 1 is updated to -1 to become state C.
[0032] When MCMC is run by introducing a parameter T that controls the above state transitions and gradually suppressing the transitions, the state asymptotically converges to the lowest energy state (ground state). Simulated annealing (hereafter referred to as SA) is a method that utilizes this to find the optimal solution to a minimization problem. The parameter T is called the temperature parameter, corresponding to real-world annealing.
[0033] When applying MCMC or SA to an interaction model, the values of variables are probabilistically updated based on equation 2 or the like. Here, variables that do not have nonlinear terms between them and are not subject to constraints are independent in the energy function equation 1 and therefore in the probability density function equation 2, and it is possible to simultaneously apply state transitions based on equation 2. Therefore, by updating independent variables in parallel, it is possible to speed up the processing of MCMC or SA.
[0034] The upper part of Figure 2 shows the interaction relationships between variables in the interaction model as a complete graph. The dotted-line enclosures indicate sets of variables with one k-hot constraint imposed. The square nodes correspond to binary variables in {-1, +1}, and the circle nodes correspond to continuous variables in [-1, +1]. The lower part shows an interaction model in which the binary variables in the interaction model shown in the upper part are relaxed to continuous values so that they can take real values.
[0035] In the upper part of Figure 2, the number of variables in the interaction model is N=6. Because there are interactions between all variables, it is a complete graph (fully connected graph). When a combinatorial optimization problem is expressed as a problem of searching the ground state of an interaction model, the interactions between variables often result in individual variables being coupled to many other variables, and therefore they are not mutually independent in the energy function as shown in the upper part of Figure 2.
[0036] Furthermore, there may be constraints on the possible states of the interaction model. A typical one is the k-hot constraint that requires k out of a set of n binary variables to be 1 (k is an integer satisfying 0 < k < n). Here, we consider the case where there is no overlap in variables among the k-hot constraints. For m = 1, ..., M k m -hot constraints (k m is an integer satisfying 0 < k m < n m ), the binary variables s im ∈ {-1, 1} belonging to each constraint are restricted as follows.
[0037]
Equation
[0038] The upper diagram in Figure 2 shows the k-hot constraint for the case of n = 4. Due to the imposition of the k-hot constraint, the possible states of a set of n binary variables are limited to n C k states, so it becomes impossible to independently transition n states so as not to transition to other states. Therefore, for example, in Patent Document 1, n C k it is proposed to improve the processing speed by parallelizing the calculation of the energy function for the possible state transitions.
[0039] As described above, due to the fact that variables are not independent in the energy function and the constraint conditions, it becomes impossible to simultaneously perform the probabilistic update of variables that satisfies the theoretical background required by MCMC. Therefore, it is difficult to speed up the processing of MCMC and SA. Especially when k ≥ 2 in the k-hot constraint, the candidates for the next state are L C k(L: dimension of the k-hot vector of interest), it is difficult to efficiently update the state while changing as many variable values as possible. Conversely, if it is possible to simultaneously update multiple spins while satisfying the constraints by constructing a set of variables that are independent in the energy function and can be moved simultaneously while satisfying the constraints, it will be possible to speed up the processing of MCMC and SA.
[0040] Patent documents 2 and 3 propose a method for efficiently finding the ground state of a fully connected interaction model when there are no constraints between variables, by solving the ground state search problem of an alternative interaction model (hereinafter referred to as an alternative model) that has a complete bipartite graph structure.
[0041] In light of this background, this embodiment describes a method for solving the ground state search problem of an interaction model having a set of binary variables to which one or more independent k-hot constraints are imposed, through the ground state search problem in an alternative model in which the interaction relationship has a complete bipartite graph structure.
[0042] First, in the interaction model, the binary variable s i We focus on the k-hot constraint imposed on ∈{-1,1}. Let S be the set of indices of the variables to which the m-th k-hot constraint is imposed. m Let i,j∈S be the matrix J in Eq. m The leading submatrix obtained by extracting the element in the i-th row and j-th column is J. m In this case, -J m The maximum eigenvalue of λ m as the constant d m =max{0,λ m}. In this case, the i-th diagonal element of J is i∈S m Then d m (i.e., J ii ←J ii +d m Let J' be the matrix where Equation 3 is satisfied and all variables are treated as continuous variables s'. i ∈[-1,1] and relaxed
[0043]
number
[0044] It can be shown that the ground state of the resulting interaction model coincides with the optimization solution of the original k-hot constrained mixed binary quadratic problem.
[0045] Generally, optimization problems that include integer variables are called integer programming problems, but when there is a mixture of integer-valued and real-valued variables, it is called a mixed integer programming problem. A mixed integer programming problem whose objective function is a quadratic function is called a mixed integer quadratic programming problem, but in this specification, a mixed integer quadratic programming problem that specifically includes a mixture of binary-valued and real-valued variables will be called a mixed binary quadratic programming problem.
[0046] Therefore, in the following, we will write the matrix J' as matrix J, the energy function H' as H, relax all binary variables to continuous variables, and define the domain of all variables as s i ∈[-1,1]. Also, let the ground state of H(s) be s * Let's say.
[0047] auxiliary variable vector v∈R N (where R is the set of real numbers)
[0048]
number
[0049] Here, I is an N-dimensional identity matrix. If the magnitude of the coupling w is greater than max{0,λ / 2} for the maximum eigenvalue λ of the interaction matrix -J, then J+2wI is a positive definite matrix, so in the ground state of H(s,v), v * = 0. In addition, the domain of s and v in Equation 5 is
[0050]
number
[0051] Even if we restrict it to , the ground state of H(s,v) is s * ,v * =0.
[0052] For s and v in Equation 6, we introduce N-dimensional variable vectors x and y defined as s = (x + y) / 2 and v = (xy) / 2.
[0053]
number
[0054] Furthermore, if we re-express Equation 5 in terms of x and y, we get the following:
[0055]
number
[0056] The upper diagram in Figure 3 is the same as the interaction model in the lower diagram in Figure 2. The lower diagram in Figure 3 is an alternative model in which the adjacency matrix between nodes has a complete bipartite graph structure, which is used to alternatively perform ground state search for the interaction model represented in the upper diagram.
[0057] The interaction relationship of Equation 8 can be expressed as an undirected graph as shown in the lower diagram of Figure 3. Also, as shown in the lower diagram of Figure 3, by imposing the same k-hot constraints on x and y as on s in Equation 8, s=(x+y) / 2 also satisfies the k-hot constraints. Therefore, the ground state of Equation 6 when the k-hot constraints are imposed on x and y is defined as x * ,y * Then s * =x * =y * holds true.
[0058] In the energy function expressed by Equation 8, each x in the variable set x(y) i (y i ) are independent of each other, as shown in the bottom part of Figure 3.
[0059] Figure 4 shows a method for updating variables in the interaction model with the complete bipartite graph structure of Figure 3. As shown in Figure 4, MCMC updates can be performed independently for each set of unconstrained variables or variables with independent k-hot constraints, so parallel processing can be used to speed up MCMC and SA processing. The following explains this parallel processing of updates.
[0060] The states are transitioned according to Gibbs sampling for the energy function of Equation 8. As the temperature parameter T, the states {x, y} of the interaction model on the right-hand side of Equation 8 are assumed to appear with a probability density function p(x, y) based on the Boltzmann distribution.
[0061]
number
[0062] Z is the normalization factor of the Boltzmann distribution, also known as the partition function. By applying Gibbs sampling to this Boltzmann distribution, it is possible to probabilistically update x and y. As an example, we will show how to update x given y. To update y given x, simply apply the formula with x and y swapped. We also define the following vector.
[0063]
number
[0064] First, the variable x that is not subject to the k-hot constraint i Here's how to update x i The probability density function of follows from Equation 9:
[0065]
number
[0066] Here, the probability density function of the truncated normal distribution with a lower limit l and an upper limit u of the random variable X is expressed by the following formula:
[0067]
number
[0068] Here, φ(·) is the probability density function of the standard normal distribution, and Φ(·) is the cumulative distribution function of the standard normal distribution function. A method for randomly sampling states from such a truncated normal distribution is described in Non-Patent Document 1. Therefore, from Equation 11, x i can be obtained by random sampling from a truncated normal distribution with the following parameters:
[0069]
number
[0070] Next, we will show how to update variables that are subject to k-hot constraints. As mentioned above, SA processing can be made more efficient by updating as many variables as possible in parallel. Therefore, here we will show a method that allows each pair to be updated independently while satisfying the k-hot constraints by creating pairs within the set of variables that are subject to k-hot constraints, as shown in Figure 4.
[0071] The solid line in Figure 5A shows the x when the sum of the pair of variables is held constant. i This is the range of motion. FIG. 5B is a graph showing the probability density function followed by the variables over the range of motion of FIG. 5A. Two variables x from the k-hot constraints i ,x j Extract these variables as x i +x j The k-hot constraint is always satisfied within the range of x, where R is kept constant. i The range of movement of x is as shown in Figure 5A. i The probability density function of follows from Equation 8:
[0072]
number
[0073] In other words, it becomes a truncated normal distribution as shown in Figure 5B. i can be obtained by random sampling from a truncated normal distribution with the following parameters:
[0074]
number
[0075] Also x j is a randomly sampled x i Using x j =Rx i Just update it as shown here. i ,x j Since the pairwise update of x does not depend on other elements of x, all pairs can be updated in parallel as shown in Figure 4.
[0076] One way to create the above pairs is to randomly create floor(n / 2) pairs (the floor function floor(·) means the largest integer less than ·) from the n variables with k-hot constraints. Another method is to create pairs so that the value of the variable changes as much as possible in the update based on Equation 14. As shown in the bottom of Figure 5, x i To avoid localization of the distribution of g i ,g j Since it is desirable that is close to 0, it is possible to create pairs that will achieve this. Furthermore, as another method based on the same idea, it is also possible to determine pairs so that the value of the variance in Equation 14 becomes large.
[0077] To implement these methods, the energy calculation execution unit 815 may determine pairs according to a rule that causes the change in the variable due to the update to be equal to or greater than a predetermined threshold. i ,g jIt is conceivable to determine pairs according to a rule such that the absolute value of σ is equal to or less than a predetermined threshold. Alternatively, it is conceivable to determine pairs according to a rule such that the variance value of a truncated normal distribution from which variables are sampled is equal to or greater than a predetermined threshold.
[0078] Next, an embodiment of an information processing device that performs the above ground state search will be described. 6 includes a processor 11, a main memory device 12, an auxiliary memory device 13, an input device 14, an output device 15, a communication device 16, one or more arithmetic units 20, and a system bus 5 that communicatively connects these devices. The information processing device 10 may be realized, for example, using virtual information processing resources such as a cloud server, part or all of which is provided by a cloud system. The information processing device 10 may also be realized, for example, by multiple information processing devices that are communicatively connected and operate in cooperation with each other.
[0079] The processor 11 is configured using, for example, a CPU (Central Processing Unit) or an MPU (Micro Processing Unit). The main storage device 12 is a device that stores programs and data, and is, for example, a ROM (Read Only Memory) (SRAM (Static Random Access Memory), NVRAM (Non Volatile RAM), Mask ROM (Mask Read Only Memory), PROM (Programmable ROM), etc.), a RAM (Random Access Memory) (DRAM (Dynamic Random Access Memory)), etc. The auxiliary storage device 13 is, for example, a hard disk drive, a flash memory, an SSD (Solid State Drive), an optical storage device (CD (Compact Disc), DVD (Digital Versatile Disc), etc.), etc. The programs and data stored in the auxiliary storage device 13 are loaded into the main storage device 12 as needed.
[0080] The input device 14 is a user interface that accepts information input from a user, and is, for example, a keyboard, a mouse, a card reader, a touch panel, etc. The output device 15 is a user interface that provides information to a user, and is, for example, a display device (LCD (Liquid Crystal Display), a graphics card, etc.) that visualizes various information, an audio output device (speaker), a printer, etc. The communication device 16 is a communication interface that communicates with other devices, and is, for example, a NIC (Network Interface Card), a wireless communication module, a USB (Universal Serial Interface) module, a serial communication module, etc.
[0081] The arithmetic device 20 is a device that executes a ground state search. The arithmetic device 20 may take the form of an expansion card attached to the information processing device 10, such as a GPU (Graphics Processing Unit). The arithmetic device 20 is configured with hardware such as a CMOS (Complementary Metal Oxide Semiconductor) circuit, an FPGA (Field Programmable Gate Array), or an ASIC (Application Specific Integrated Circuit). The arithmetic device 20 includes a control device, a storage device, an interface for connecting to the system bus 5, and transmits and receives commands and information to and from the processor 11 via the system bus 5. The arithmetic device 20 may be communicably connected to another arithmetic device 20 via a communication line, for example, and may operate in cooperation with the other arithmetic device 20. The functions realized by the arithmetic device 20 may be realized, for example, by causing a processor (CPU, GPU, etc.) to execute a program.
[0082] 7 is a diagram illustrating the operating principle of the arithmetic device 20, and is a block diagram of a circuit (hereinafter referred to as an arithmetic circuit 700) that constitutes the arithmetic device 20. The arithmetic circuit 700 realizes a function equivalent to the update formula. The operating principle of the arithmetic device 20 will be described below with reference to this diagram.
[0083] As shown in the figure, the arithmetic circuit 700 includes a nonlinear coefficient memory 711, a linear coefficient memory 712, a D-th variable set memory 713.D (D=1, 2), a D-th variable set pair list memory 714.D (D=1, 2), a product-sum arithmetic unit 715, and a mathematical function arithmetic unit 716.
[0084] The nonlinear coefficient memory 711 stores information representing a real symmetric matrix J consisting of nonlinear coefficients that represent interactions between nodes and information representing couplings w (see Equation 5). As mentioned above, the interaction matrix is a real symmetric matrix, and therefore, this symmetry can be used to reduce the amount of nonlinear coefficient memory 711 used. The linear coefficient memory 712 stores information representing the linear coefficient vector h (see Equation 5). The D-th variable set memory 713.D (D=1, 2) stores information on N-dimensional vectors indicating the state of the D-th variable set of the substitution model described above (see formula 5). The D-th variable set memory 713.D is configured to be able to store real values (continuous values).
[0085] Signals SW, SR, and ST are input to the arithmetic circuit 700. A mathematical function arithmetic unit 716 outputs a signal SP. The signal SW is a signal that periodically repeats the integers 1 and 2, and designates the first variable set memory and the second variable set memory. The signal SR is a signal representing a vector in which each element is a mutually independent random number. The random numbers used in the stochastic update described in Non-Patent Document 1 are input. The signal ST inputs the temperature parameter T at SA.
[0086] As described above in the explanation of Equation 5, the coupling w is set based on the eigenvalue information of the matrix J. The numerical evaluation may be performed outside the arithmetic unit 20, for example, by the processor 11. Alternatively, the calculation may be performed within the arithmetic unit 20. For example, when calculating the maximum eigenvalue by the power method or the like, this involves repeatedly executing a matrix-vector product, and the product-sum calculation unit 715 can be used.
[0087] A nonlinear coefficient memory 711, a linear coefficient memory 712, a D-th variable group memory 713.D, and a signal SW are input to a product-sum calculation unit 715. Equation 10 includes product-sum calculations of a matrix and a vector, and the calculation is executed and output under the corresponding relationship of the variables specified by the signal SW.
[0088] The random sampling of the truncated normal distribution in Non-Patent Document 1 includes the calculation of mathematical functions such as the output of a product-sum operation and an error function for the uniform random number r, and this part is executed by the mathematical function calculation unit 716. The signal SW specifies which equation corresponds to which calculation to execute. The uniform random number r and the temperature parameter T are input via signals SR and ST, respectively. The output values, which are the values of the variables updated according to MCMC, are output via signal SP.
[0089] Although the scope of disclosures in Patent Documents 2 and 3 will not be described in detail in this specification, state transitions of variables are sequentially executed through stochastic state transitions according to a simulated annealing algorithm. Here, calculations of update formulas can be executed independently for independent variables or pairs in each variable group. In other words, the next state is determined stochastically by MCMC so that variables without k-hot constraints are independently determined, and variables with k-hot constraints are determined independently so that the sum of the variable sets is constant. This allows variable updates to be executed in parallel by multiple arithmetic circuits 700, thereby realizing high-speed MCMC. Furthermore, even when k-hot constraints are imposed, the degree of freedom in variable updates is large.
[0090] 8 shows the main functions (software configuration) of the information processing device 10. As shown in the figure, the information processing device 10 includes a storage unit 800, a model coefficient setting unit 811, a coupling strength calculation unit 812, a variable value initialization unit 813, a temperature parameter control unit 814, an energy calculation execution unit 815, and a variable value reading unit 816. These functions are realized by the processor 11 reading and executing a program stored in the main storage device 12, or by hardware included in the calculation device 20. In addition to the above functions, the information processing device 10 may also include other functions, such as an operating system, a file system, a device driver, a DBMS (DataBase Management System), etc.
[0091] Of the above functions, the storage unit 800 stores MBQP-format problem data 801 and a calculation device control program 802 in the main storage device 12 or the auxiliary storage device 13. The MBQP-format problem data 801 is data in which a combinatorial optimization problem is input in a predetermined description format, and holds constraint conditions such as nonlinear coefficients, linear coefficients, and k-hot constraints. The MBQP-format problem data 801 is set, for example, by a user via a user interface (input device, output device, communication device, etc.).
[0092] The calculation device control program 802 is a program that is executed when the energy calculation execution unit 815 controls the calculation device 20, or that is loaded by the energy calculation execution unit 815 into each calculation device 20 and executed by the calculation device 20.
[0093] When the MBQP-style problem data 801 is a mixed binary quadratic problem, the energy calculation execution unit 815 continuously relaxes the binary variables of the MBQP-style problem data 801 with k-hot constraints based on the calculation device control program 802. That is, the variables are made to take on real values. Furthermore, variable pairs are created within the k-hot vector, and control is performed so that state updates are performed while keeping their sum constant. The storage unit 800 may previously store MBQP-style problem data 801 that has been continuously relaxed.
[0094] After obtaining the optimal solution, the energy calculation execution unit 815 performs control to convert the continuous variables back to binary variables. Conversion between binary variables and continuous variables can be performed based on the disclosure of Patent Document 3.
[0095] The model coefficient setting unit 811 sets the non-linear coefficient memory 711 and the linear coefficient memory 712 based on the MBQP format problem data 801 .
[0096] The connection strength calculation unit 812 sets the value of the connection w based on the information in the nonlinear coefficient memory 711 and the eigenvalues of the interaction matrix.
[0097] The variable value initialization unit 813 initializes values stored in the variable memory of the arithmetic unit 20. For example, the values are set by uniform random numbers between -1 and +1 while satisfying each k-hot constraint.
[0098] The temperature parameter control unit 814 controls the temperature parameter T in the SA.
[0099] The energy calculation execution unit 815 performs a ground state search (hereinafter referred to as energy calculation) for the interaction model according to SA. For this purpose, the energy calculation execution unit 815 also performs overall control to cause other units, such as the model coefficient setting unit 811, the coupling strength calculation unit 812, the variable value initialization unit 813, the temperature parameter control unit 814, and the variable value reading unit 816, to execute necessary processes.
[0100] When SA is executed by the energy calculation execution unit 815, the variable value reading unit 816 reads the value stored in the variable memory and outputs the read value to the output device 15 or the communication device 16, thereby completing the ground state search.
[0101] 9 is a flowchart for explaining the processing performed by the information processing device 10 when searching for the ground state of an interaction model (hereinafter referred to as ground state search processing S900). The ground state search processing S900 will be explained below with reference to the same figure. Note that, in the following, the letter "S" added before the reference numeral means a processing step. The ground state search processing S900 is started, for example, by receiving an instruction from the user via the input device 14.
[0102] First, the model coefficient setting unit 811 sets values in the nonlinear coefficient memory 711 and the linear coefficient memory 712 of the arithmetic device 20 (arithmetic circuit 700) based on the MBQP format problem data 801 (S911). The memory values can also be set or edited by the user via a user interface (realized by, for example, the input device 14, the output device 15, the communication device 16, etc.).
[0103] Next, the connection strength calculation unit 812 sets a correction value d based on the k-hot constraint and the matrix J stored in the nonlinear coefficient memory 711, and stores the correction value d in the nonlinear coefficient memory 711 again. As described above, this calculation may be performed in the arithmetic unit 20 or the processor 11 (S912). The connection strength calculation unit 812 uses the constant d as the correction value d, as described above. m =max{0,λ m}.
[0104] Next, the coupling strength calculation unit 812 calculates the matrix J'(J ii ←J ii +d m The coupling w is set based on the matrix J′) and is newly stored in the nonlinear coefficient memory 711. As described above, this calculation may be performed in the arithmetic unit 20 or the processor 11 (S913).
[0105] Next, the variable value initialization unit 813 initializes the values stored in the variable memory (S914). Since the energy calculation execution unit 815 redefines the binary variables defined in the MBQP format problem data 801 as continuous variables, the variable values are continuous.
[0106] Next, the energy calculation execution unit 815 updates the pair list, which is a list of variables to be updated in pairs based on random or probability density function information, and stores it in the D-th variable set pair list memory 714. The calculation circuit 700 updates the value of the D-th variable set memory 713 based on the result of the energy calculation execution unit 815. This operation is applied alternately to the first variable set and the second variable set (S915).
[0107] The method by which the arithmetic device 20 (arithmetic circuit 700) performs state updates (variable updates) and transitions to the ground state may follow, for example, the methods of Patent Document 2 and Patent Document 3. However, in this embodiment, after continuously relaxing the variables, control is performed to keep the sum of the values of a predetermined variable pair constant based on the contents of the Dth variable group pair list memory 714. Since variable pairs can be updated independently without affecting other variables, all pairs can be updated in parallel.
[0108] Next, the energy calculation execution unit 815 determines whether the SA termination condition is met (for example, whether the state update has been executed while changing the temperature parameter T a predetermined number of times) (S916). If the energy calculation execution unit 815 determines that the SA termination condition is met (S916: YES), the process proceeds to S917. On the other hand, if the energy calculation execution unit 815 determines that the stop condition is not met (S916: NO), the process returns to S915.
[0109] Next, the variable value reading unit 917 reads the value stored in the variable memory and stores it as the result of the ground state search (S917), and the ground state search process S900 ends.
[0110] Based on the results of the ground state search, the energy calculation execution unit 815 performs processing to convert the continuous variables of the solution back into binary variables. For example, a binary variable of -1 or +1 may be obtained based on the sign of the continuous variable.
[0111] 10 is a block diagram showing a detailed configuration example of the arithmetic unit 20, and is a block diagram showing a circuit configuration example in which SRAM technology is applied to this embodiment to form a semiconductor integrated circuit. A unit that executes the operation of the arithmetic circuit 700 constitutes an array unit 1001. Such a configuration can be manufactured by applying semiconductor manufacturing technology.
[0112] An example configuration will be described in Fig. 10. Data stored in the nonlinear coefficient memory 711 and the linear coefficient memory 712 is set by a model coefficient setting unit 911. The matrix J is stored in the nonlinear coefficient memory 711, and the linear coefficient vector h is stored in the linear coefficient memory 712, but these are used in common by all units to reduce the circuit size. The signal lines for this purpose are omitted in Fig. 10. In principle, each unit may be provided with its own memory.
[0113] The SRAM interface 1002 writes to and reads from the variable set memory. After the processing in the arithmetic circuit 700 is completed, the values of the variables read are sent to the variable value reading unit 816. The variable value reading unit 816 outputs the results of the ground state search by storing and outputting the read values as appropriate.
[0114] The controller 1003 initializes each memory and reports the completion of the calculation process of each unit in response to instructions from the energy calculation execution unit 815 .
[0115] Although one embodiment has been described above in detail, it goes without saying that the present invention is not limited to the above embodiment and can be modified in various ways without departing from the spirit of the present invention. For example, the above embodiment has been described in detail to clearly explain the present invention, and the present invention is not necessarily limited to those having all of the described configurations. Furthermore, some of the configurations of the above embodiment can be added, deleted, or replaced with other configurations.
[0116] Furthermore, the above-described configurations, functional units, processing units, processing means, etc. may be partially or entirely implemented in hardware, for example, by designing them as integrated circuits. The above-described configurations, functions, etc. may also be implemented in software, with a processor interpreting and executing a program that implements each function. Information such as the programs, tables, and files that implement each function can be stored in a memory, a hard disk, a recording device such as an SSD (Solid State Drive), an IC card, an SD card, a DVD, or other recording media.
[0117] In addition, in the above figures, the control lines and information lines shown are those that are considered necessary for explanation, and do not necessarily show all the control lines and information lines that are actually implemented. For example, it may be considered that almost all components are actually connected to each other.
[0118] The above-described layout of the various functional units, processing units, and databases of the information processing device 10 is merely an example. The layout of the various functional units, processing units, and databases can be changed to an optimal layout in terms of the performance, processing efficiency, communication efficiency, etc. of the hardware and software included in the information processing device 10.
[0119] Furthermore, the configuration (schema, etc.) of the database that stores the various types of data described above can be flexibly changed from the viewpoint of efficient use of resources, improved processing efficiency, improved access efficiency, improved search efficiency, and the like.
[0120] As described above, using an information processing device, a mixed binary quadratic optimization problem with k-hot constraints for binary variables is reduced to a quadratic programming problem with k-hot constraints by adding correction values to the nonlinear coefficients based on the k-hot constraints and information on the nonlinear coefficients. Furthermore, the optimal solution search for this quadratic programming problem is performed by ground-state search for an interaction model with a complete bipartite graph structure when the nonlinear coefficients are expressed as an adjacency matrix. This ground-state search is performed by updating variables according to simulated annealing based on the energy function of the interaction model and information on the variables. This method improves processing speed by solving a mixed binary quadratic optimization problem with k-hot constraints for binary variables using continuous-value relaxation.
[0121] As described above, we can provide a method for efficiently searching for the optimal solution to a mixed binary quadratic optimization problem that has one or more independent k-hot constraints on binary variables, based on continuous value relaxation of the binary variables and parallel processing using this. [Explanation of symbols]
[0122] 5 system bus, 10 information processing device, 11 processor, 12 main memory, 13 auxiliary memory, 14 input device, 15 output device, 16 communication device, 12 main memory, 20 arithmetic unit, 711 nonlinear coefficient memory, 712 linear coefficient memory, 713.D(D=1,2) Dth variable group memory, 714.D(D=1,2) Dth variable group pair list memory, 715 product-sum calculation unit, 716 mathematical function calculation unit, 800 memory unit, 801 MBQP format problem data, 802 arithmetic unit control program, 811 model coefficient setting unit, 812 coupling strength calculation unit, 813 variable value initialization unit, 814 temperature parameter control unit, 815 energy calculation execution unit, 816 variable value reading unit, 1001 unit array, 1102 SRAM interface, 1103 controller
Claims
1. An optimization method for solving a quadratic programming problem in which one or more independent k-hot constraints (k is a positive natural number) are imposed on binary variables, using an information processing device including a processor, a storage device, an input device, and an output device, comprising: the information processing device adds a correction value to a nonlinear coefficient of the binary variable to which the k-hot constraint is imposed, thereby relaxing the binary variable to a continuous value; The information processing device performs a state transition so that the sum of the set of continuous-valued variables to which the k-hot constraint is imposed becomes constant, thereby performing a solution search while satisfying the k-hot constraint. Optimization methods.
2. The correction value of the nonlinear coefficient of the binary variable to which the k-hot constraint is imposed is The information processing device determines nonlinear coefficients between binary variables to which a k-hot constraint is imposed based on eigenvalues of the extracted leading submatrix. The optimization method of claim 1 .
3. The information processing device includes: Two variable sets x and y each having N variables are stored, storing an N-th order real symmetric matrix J defined by nonlinear coefficients between variables of the quadratic programming problem and a vector h defined by linear coefficients; The i-th pair of variables x of the two variable sets i , y i During this time, a coupling strength w is calculated based on information on the eigenvalues of the N-th order real symmetric matrix J. The optimization method of claim 1 .
4. The solution search is The N-th order real symmetric matrix J acts as an adjacency matrix between the two sets of variables x and y, the vector h acts as a bias coefficient for x and y, and a ground state search for an interaction model is performed such that an undirected graph expressing x and y as nodes has a complete bipartite graph structure.
4. The optimization method of claim 3.
5. The ground state search This is done by sequentially executing the state transitions of variables through stochastic state transitions according to the simulated annealing algorithm.
5. The optimization method of claim 4.
6. The state transition is is performed simultaneously on a plurality of variables belonging to the set of variables x or y, A variable without k-hot constraints is The variables to which the k-hot constraint is imposed are set so that the sum is constant. The next state is determined probabilistically using the Markov Chain Monte Carlo method.
6. The optimization method of claim 5.
7. The quadratic programming problem is a mixed binary quadratic programming problem. The optimization method of claim 1 .
8. the information processing device randomly determines the set of continuous variables; The optimization method of claim 1 .
9. the information processing device determines the set of continuous-value variables so that a change in the variable in the state transition is equal to or greater than a predetermined threshold; The optimization method of claim 1 .
10. An information processing device comprising a processor, a storage device, an input device, an output device, and an arithmetic unit, for solving a quadratic programming problem in which one or more independent k-hot constraints (k is a positive natural number) are imposed on binary variables, an energy calculation execution unit that performs continuous value relaxation on the binary variables to which the k-hot constraint is imposed; a connection strength calculation unit that calculates a correction value to be added to a nonlinear coefficient of the binary variable to which the k-hot constraint is imposed, Under the control of the energy calculation execution unit, the calculation device performs a state transition so that the sum of the set of continuous-valued variables to which the k-hot constraint is imposed becomes constant, thereby performing a solution search while satisfying the k-hot constraint. Information processing device.
11. The connection strength calculation unit calculating a correction value of the nonlinear coefficient of the binary variables to which the k-hot constraint is imposed based on the eigenvalue of the leading submatrix that extracts the nonlinear coefficient between the binary variables to which the k-hot constraint is imposed; The information processing device according to claim 10.
12. The arithmetic unit is configured by a semiconductor integrated circuit, and includes a variable memory for storing two sets of variables x and y, each having N variables, and a nonlinear coefficient memory for storing an N-th order real symmetric matrix J defined by nonlinear coefficients between the variables of the quadratic programming problem, The connection strength calculation unit calculates the i-th variable pair x of the two variable sets based on the information in the nonlinear coefficient memory. i , y i During this time, a coupling strength w is calculated based on information on the eigenvalues of the N-th order real symmetric matrix J. The information processing device according to claim 11.
13. the calculation device includes a linear coefficient memory for storing a vector h that is a bias coefficient for x and y; The N-th order real symmetric matrix J is applied between the two sets of variables x and y as an adjacent matrix, Apply the vector h to x and y, Execute a ground state search for an interaction model such that an undirected graph representing x and y as nodes has a complete bipartite graph structure.
13. The information processing device according to claim 12.
14. the computing device performs the ground state search by sequentially executing state transitions of variables through stochastic state transitions in accordance with a simulated annealing algorithm; 14. The information processing device according to claim 13.
15. In the state transition, is performed simultaneously on a plurality of variables belonging to the set of variables x or y, A variable without k-hot constraints is The variables to which the k-hot constraint is imposed are set so that the sum is constant. Probabilistic determination of the next state using the Markov Chain Monte Carlo method 15. The information processing device according to claim 14.
Citation Information
Patent Citations
Optimization device and method for controlling optimization device
JP2020064535A
Information processing system
JP2022083776A
Information processing device, calculation device, and information processing method
WO2019216277A1
Systems and methods of hybrid algorithms for solving discrete quadratic models
WO2021126773A1
Computation system, information processing device, and optimum solution search processing method
WO2021220445A1