Information processing device, information processing system, information processing method, storage medium and computer program product

Abstract

Provided are an information processing device, an information processing system, an information processing method, a storage medium, and a computer program product for calculating a solution to a combinatorial optimization problem within a practical time. An information processing device as an embodiment of the present invention includes a storage unit and a processing circuit. The storage unit is configured to store a first variable as an element of a first vector and a second variable as an element of a second vector. The processing circuit is configured to update the first variable based on the corresponding second variable, weight the first variable with a first coefficient and add the corresponding second variable, calculate a problem item using a plurality of the first variables, add the second variable to the problem item, calculate a first correction item including the product of a constraint item and a second coefficient, add the first correction item to the second variable, and increase the absolute values ​​of the first coefficient and the second coefficient according to the number of updates. The constraint item is based on a constraint function representing a constraint condition and has the first variable as an independent variable.

Description

Technical Field

[0001] Embodiments of the present invention relate to an information processing device, an information processing system, an information processing method, a storage medium, and a computer program product. Background Art

[0002] The combinatorial optimization problem is a problem of selecting the combination that best suits the purpose from among multiple combinations. Mathematically, the combinatorial optimization problem is reduced to the problem of maximizing or minimizing a function with multiple discrete variables called the "objective function". The combinatorial optimization problem is a universal problem in various fields such as finance, logistics, transportation, design, manufacturing, and life sciences, but the optimal solution is not necessarily found due to the so-called "combinatorial explosion" in which the number of combinations increases at the order of an exponential function of the problem size. In addition, there are many cases where it is difficult to even obtain an approximate solution that is close to the optimal solution.

[0003] In order to solve problems in various fields and promote social innovation and scientific and technological progress, it is required to develop technology that can calculate solutions to combinatorial optimization problems in a practical time.

[0004] Prior art literature

[0005] Patent Literature

[0006] Patent Document 1: Japanese Patent Application Publication No. 2017-73106

[0007] Non-patent literature

[0008] Non-patent literature 1: H. Goto, K. Tatsumura, ARDixon, Sci. Adv. 5, eaav2372 (2019)

[0009] Non-patent document 2: H. Goto, Sci. Rep. 6, 21686 (2016)

[0010] Non-Patent Literature 3: Tsuchiya, Nishiyama, and Tsujita: Approximate solution to combinatorial optimization problems using branching characteristics

[0011] Non-Patent Literature 4: Tsuchiya, Nishiyama, and Tsujita: Analysis of the Deterministic Annealing Algorithm Summary of the invention

[0012] Technical problem to be solved by the invention

[0013] The object of the present invention is to provide an information processing device, an information processing system, an information processing method, a storage medium and a computer program product for calculating the solution of a combinatorial optimization problem within a practical time.

[0014] Means for solving technical problems

[0015] An information processing device as an embodiment of the present invention includes a storage unit and a processing circuit. The storage unit is configured to store a first variable as an element of a first vector and a second variable as an element of a second vector. The processing circuit is configured to update the first variable based on the corresponding second variable, weight the first variable with a first coefficient and add the corresponding second variable, calculate a problem item using a plurality of the first variables, add the second variable to the problem item, calculate a first correction item including the product of a constraint item and a second coefficient, add the second variable to the first correction item, and increase the absolute values ​​of the first coefficient and the second coefficient according to the number of updates. The constraint item is based on a constraint function and has the first variable as an independent variable. BRIEF DESCRIPTION OF THE DRAWINGS

[0016] Figure 1 It is a diagram showing a configuration example of an information processing system.

[0017] Figure 2 This is a block diagram showing an example configuration of a management server.

[0018] Figure 3 This is a diagram showing an example of data stored in the storage unit of the management server.

[0019] Figure 4 This is a block diagram showing an example of the configuration of a calculation server.

[0020] Figure 5 The diagram shows an example of data stored in a storage device of a computing server.

[0021] Figure 6 This is a flowchart showing an example of processing when a solution of a simulation branching algorithm is calculated by time advancement.

[0022] Figure 7 is a table showing examples of equality constraints.

[0023] Figure 8 is a table showing examples of inequality constraints.

[0024] Fig. 9 is a table showing examples of inequality constraints.

[0025] Fig.10 This is a flowchart showing an example of a process of solving an extended Hamiltonian including a first correction term.

[0026] Fig.11 : is a flowchart showing an example of a solution process of an extended Hamiltonian including a second correction term.

[0027] Fig.12is the coefficient λ m Flowchart of an example of a solution process in a case where a part of the update process is skipped.

[0028] Fig.13 This is a diagram schematically showing an example of a multi-processor configuration.

[0029] Fig.14 This is a diagram schematically showing an example of a configuration using a GPU.

[0030] Fig.15 is a flow chart showing an example of the overall processing performed to solve a combinatorial optimization problem. DETAILED DESCRIPTION

[0031] Below, while referring to the attached Figure 1 In the drawings, the same components are denoted by the same reference numerals, and description thereof will be appropriately omitted.

[0032] Figure 1 It is a block diagram showing a configuration example of the information processing system 100 . Figure 1 The information processing system 100 includes a management server 1, a network 2, computing servers (information processing devices) 3a to 3c, cables 4a to 4c, a switch 5, and a storage device 7. Figure 1 , represents a client terminal 6 capable of communicating with the information processing system 100. The management server 1, the computing servers 3a to 3c, the client terminal 6, and the storage device 7 can communicate data with each other via the network 2. For example, the computing servers 3a to 3c can save data to the storage device 7 or read data from the storage device 7. The network 2 is, for example, the Internet formed by interconnecting multiple computer networks. The network 2 can use wired, wireless, or a combination of these as a communication medium. In addition, as an example of a communication protocol used in the network 2, there is TCP / IP, but there is no limitation on the type of communication protocol.

[0033] In addition, the computing servers 3a~3c are connected to the switch 5 via cables 4a~4c, respectively. The cables 4a~4c and the switch 5 form an intermediate connection between the computing servers. The computing servers 3a~3c can communicate data with each other via the intermediate connection. The switch 5 is, for example, an Infiniband switch. The cables 4a~4c are, for example, Infiniband cables. However, a wired LAN switch / cable can be used instead of an Infiniband switch / cable. There is no particular limitation on the communication standard and communication protocol used in the cables 4a~4c and the switch 5. Examples of the client terminal 6 include a notebook PC, a desktop PC, a smart phone, a tablet computer, a car terminal, and the like.

[0034] In solving the combinatorial optimization problem, the processing can be parallelized and / or decentralized. Therefore, the computing servers 3a~3c and / or the processors of the computing servers 3a~3c can share the steps of executing a part of the computing processing, and can also perform the same computing processing in parallel for different variables. The management server 1, for example, transforms the combinatorial optimization problem input by the user into a form that can be processed by each computing server, and controls the computing server. In addition, the management server 1 obtains the calculation results from each computing server and transforms the aggregated calculation results into the solution of the combinatorial optimization problem. In this way, the user can obtain the solution to the combinatorial optimization problem. It is assumed that the solution to the combinatorial optimization problem includes an optimal solution and an approximate solution close to the optimal solution.

[0035] exist Figure 1 3 computing servers are shown in . However, the number of computing servers included in the information processing system is not limited. In addition, there is no special limitation on the number of computing servers used in solving the combinatorial optimization problem. For example, the information processing system may include one computing server. In addition, any of the multiple computing servers included in the information processing system can be used to solve the combinatorial optimization problem. In addition, the information processing system may include more than hundreds of computing servers. The computing server can be a server located in a data center or a desktop PC located in an office. In addition, the computing server can also be multiple types of computers located in different locations. There is no special limitation on the type of information processing device used as a computing server. For example, the computing server can be a general-purpose computer, a dedicated electronic circuit, or a combination thereof.

[0036] Figure 2 It is a block diagram showing a configuration example of the management server 1 . Figure 2 The management server 1 is, for example, a computer including a central processing unit (CPU) and a memory. The management server 1 includes a processor 10, a storage unit 14, a communication circuit 15, an input circuit 16, and an output circuit 17. The processor 10, the storage unit 14, the communication circuit 15, the input circuit 16, and the output circuit 17 are connected to each other via a bus 20. The processor 10 includes a management unit 11, a conversion unit 12, and a control unit 13 as internal components.

[0037] The processor 10 is an electronic circuit that performs calculations and controls the management server 1. As the processor 10, for example, a CPU, a microprocessor, an ASIC, an FPGA, a PLD or a combination thereof can be used. The management unit 11 provides an interface for operating the management server 1 via the user's client terminal 6. Examples of the interface provided by the management unit 11 include API, CLI, web pages, etc. For example, the user can input information about the combinatorial optimization problem through the management unit 11, or browse and / or download the calculated solution to the combinatorial optimization problem. The conversion unit 12 converts the combinatorial optimization problem into a form that can be processed by each computing server. The control unit 13 sends a control instruction to each computing server. After the control unit 13 obtains the calculation results from each computing server, the conversion unit 12 aggregates the multiple calculation results and converts them into a solution to the combinatorial optimization problem. In addition, the control unit 13 can also specify the processing content to be performed by each computing server or the processor in each server.

[0038] The storage unit 14 stores various data including the program of the management server 1, the data required for the execution of the program, and the data generated by the program. Here, it is assumed that the program includes both the OS and the application. The storage unit 14 can also be a volatile memory, a non-volatile memory, or a combination thereof. Examples of volatile memories include DRAM, SRAM, etc. Examples of non-volatile memories include NAND flash memory, NOR flash memory, ReRAM, or MRAM. In addition, a hard disk, an optical disk, a magnetic tape, or an external storage device can also be used as the storage unit 14.

[0039] The communication circuit 15 transmits and receives data to and from each device connected to the network 2. The communication circuit 15 is, for example, a NIC (Network Interface Card) of a wired LAN. However, the communication circuit 15 may also be another type of communication circuit such as a wireless LAN. The input circuit 16 implements data input to the management server 1. It is assumed that the input circuit 16 has, for example, a USB, PCI-Express, etc. as an external port. Figure 2 In the example, the operating device 18 is connected to the input circuit 16. The operating device 18 is a device for inputting information to the management server 1. The operating device 18 is, for example, a keyboard, a mouse, a touch panel, a voice recognition device, etc., but is not limited thereto. The output circuit 17 realizes data output from the management server 1. It is assumed that the output circuit 17 has HDMI, DisplayPort, etc. as an external port. Figure 2 In the example of , the display device 19 is connected to the output circuit 17. Examples of the display device 19 include an LCD (Liquid Crystal Display), an organic EL (Organic Electroluminescence) display, or a projector, but are not limited thereto.

[0040] The administrator of the management server 1 can use the operation device 18 and the display device 19 to perform maintenance on the management server 1. In addition, the operation device 18 and the display device 19 may be incorporated into the management server 1. In addition, the operation device 18 and the display device 19 may not necessarily be connected to the management server 1. For example, the administrator may perform maintenance on the management server 1 using a client terminal that can communicate with the network 2.

[0041] Figure 3 2 shows an example of data stored in the storage unit 14 of the management server 1. Figure 3 In the storage unit 14, problem data 14A, calculation data 14B, management program 14C, transformation program 14D, and control program 14E are stored. For example, problem data 14A includes data of a combinatorial optimization problem. For example, calculation data 14B includes calculation results collected from each calculation server. For example, management program 14C is a program that implements the functions of the management unit 11 described above. For example, transformation program 14D is a program that implements the functions of the transformation unit 12 described above. For example, control program 14E is a program that implements the functions of the control unit 13 described above.

[0042] Figure 4 FIG. 1 is a block diagram showing an example of the configuration of a computing server. Figure 4 In the figure, the structure of the computing server 3a is shown by way of example. The other computing servers may have the same structure as the computing server 3a or may have a different structure from the computing server 3a. The computing server 3a is, for example, an information processing device that performs the calculation of the first vector and the second vector alone or in sharing with other computing servers. In addition, at least one computing server may calculate the problem item between the elements of the first vector. Here, the problem item may also be calculated based on the Ising model described later. In addition, the problem item may also include the multi-body interaction described later.

[0043] For example, the first vector is based on the variable x i (i=1, 2, ..., N) as the elements. For example, the second vector is a vector with variable y i (i=1, 2, ..., N) is a vector of elements.

[0044] The calculation server 3a includes, for example, a communication circuit 31, a shared memory 32, processors 33A to 33D, a storage device 34, and a host bus adapter 35. The communication circuit 31, the shared memory 32, the processors 33A to 33D, the storage device 34, and the host bus adapter 35 are connected to each other via a bus 36.

[0045] The communication circuit 31 transmits and receives data with each device connected to the network 2. The communication circuit 31 is, for example, a NIC (Network Interface Card) of a wired LAN. However, the communication circuit 31 may also be other types of communication circuits such as a wireless LAN. The shared memory 32 is a memory accessible from the processors 33A to 33D. As an example of the shared memory 32, a volatile memory such as a DRAM and an SRAM may be cited. However, as the shared memory 32, other types of memories such as a non-volatile memory may also be used. The shared memory 32 may also be configured to store the elements of the first vector and the elements of the second vector, for example. Here, the shared memory 32 and the storage device 34 described later are examples of the storage unit of the information processing device. The processors 33A to 33D can share data via the shared memory 32. In addition, all the memories of the computing server 3a may not necessarily be configured as shared memories. For example, a part of the memory of the computing server 3a may also be configured as a local memory that can be accessed only from a certain processor.

[0046] Processors 33A to 33D are electronic circuits that perform computational processing. The processor may be, for example, a CPU (Central Processing Unit), a GPU (Graphics Processing Unit), an FPGA (Field-Programmable Gate Array), or an ASIC (Application Specific Integrated Circuit), or a combination thereof. In addition, the processor may also be a CPU core or a CPU thread. When the processor is a CPU, there is no particular limitation on the number of slots that the computing server 3a has. In addition, the processor may also be connected to other components of the computing server 3a via a bus such as PCI express.

[0047] exist Figure 4 In the example, the computing server has four processors. However, the number of processors in one computing server may be different. For example, the number and / or type of processors installed in the computing server may be different.

[0048] The information processing device is configured, for example, to i The first vector with (i=1, 2, ..., N) as its elements and the second variable y corresponding to the first variable i The second vector having (i=1, 2, ..., N) as elements is repeatedly updated. The storage unit of the information processing device may be configured to store the first variable as an element of the first vector and the second variable as an element of the second vector.

[0049] The processing circuit of the information processing device, for example, has a processing circuit configured as follows: updating the first variable based on the corresponding second variable, weighting the first variable with the first coefficient and adding it to the corresponding second variable, using multiple first variables to calculate the problem item, adding the problem item to the second variable, calculating the first correction item including the product of the constraint item and the second coefficient, adding the first correction item to the second variable, and increasing the absolute values ​​of the first coefficient and the second coefficient according to the number of updates. Here, the constraint item is based on a constraint function representing the constraint condition and has the first variable as an independent variable (Japanese: argument). Here, the above-mentioned processors 33A to 33D are an example of a processing circuit. In this way, the information processing device may also have multiple processing circuits. In this case, each processing circuit may also be configured to update at least a portion of the first vector and at least a portion of the second vector in parallel.

[0050] The following is an example in which the first coefficient and the second coefficient are positive values ​​and the values ​​of the first coefficient and the second coefficient increase corresponding to the number of updates. However, it is also possible to use negative first coefficients and second coefficients by reversing the sign of the algorithm suggested below. In this case, the values ​​of the first coefficient and the second coefficient may decrease corresponding to the number of updates. However, in either case, it can be said that the absolute values ​​of the first coefficient and the second coefficient increase corresponding to the number of updates. In addition, the problem item for processing circuit calculations may also be based on the Ising model. In addition, the problem item for processing circuit calculations may also include multi-body interactions.

[0051] Here, the case where the processing content is allocated to the processor unit is used as an example for explanation. However, the unit of computing resources for allocating the processing content is not limited. For example, the processing content can also be allocated to the computer unit.

[0052] Below, refer again Figure 4 , describing the components of a computing server.

[0053] The storage device 34 stores various data including the program of the computing server 3a, the data required for the execution of the program, and the data generated by the program. Here, it is assumed that the program includes both the OS and the application. The storage device 34 can also be configured to store the elements of the first vector and the elements of the second vector, for example. The storage device 34 can also be a volatile memory, a non-volatile memory, or a combination thereof. Examples of volatile memories include DRAM or SRAM. Examples of non-volatile memories include NAND flash memory, NOR flash memory, ReRAM, or MRAM. In addition, a hard disk, an optical disk, a magnetic tape, or an external storage device can also be used as the storage device 34.

[0054] The host bus adapter 35 implements data communication between computing servers. The host bus adapter 35 is connected to the switch 5 via the cable 4a. The host bus adapter 35 is, for example, an HCA (Host Channel Adaptor). By forming an intermediate connection that can achieve high communication throughput using the host bus adapter 35, the cable 4a, and the switch 5, the speed of parallel computing processing can be increased.

[0055] Figure 5 An example of data stored in a storage device on a computing server. Figure 5 The storage device 34 stores calculation data 34A, calculation program 34B, and control program 34C. The calculation data 34A includes data or calculation results in the process of calculation of the calculation server 3a. In addition, at least a part of the calculation data 34A can also be stored in different storage levels such as the shared memory 32, the processor cache, or the processor register. The calculation program 34B is a program that implements the calculation processing in each processor and the storage processing of the data in the shared memory 32 and the storage device 34 based on a prescribed algorithm. The control program 34C is a program that controls the calculation server 3a based on the instructions sent from the control unit 13 of the management server 1 and sends the calculation results of the calculation server 3a to the management server 1.

[0056] Next, the technology associated with solving the combinatorial optimization problem is described. As an example of an information processing device used to solve the combinatorial optimization problem, an Ising machine can be cited. An Ising machine refers to an information processing device that calculates the energy of the basic state of the Ising model. So far, the Ising model has been mainly used as a model of ferromagnetic bodies and phase transition phenomena. However, in recent years, the use of the Ising model as a model for solving combinatorial optimization problems has increased. The following formula (1) represents the energy of the Ising model.

[0057] [Formula 1]

[0058]

[0059] Here, s i 、s j is the spin, which is a binary variable that takes a value of +1 or -1. N is the spin number. h i is the local magnetic field acting on each spin. J is the matrix of coupling coefficients between spins. The matrix J is a real symmetric matrix with zero diagonal components. Therefore, J ij represents the element of row i and column j of matrix J. The Ising model in equation (1) is a quadratic equation for spin, but as described later, an extended Ising model (Ising model with multi-body interaction) including a term of degree 3 or higher for spin may also be used.

[0060] If the Ising model of formula (1) is used, the energy E Ising As the objective function, calculate the energy E Ising The solution of the Ising model is the vector of spins (s1, s2, ..., s N ) form. In particular, the energy E Ising is the vector of minimum values ​​(s1, s2, ..., s N ) is called the optimal solution. However, the solution of the Ising model calculated may not necessarily be the strict optimal solution. In the following, the Ising model will be used to calculate the energy E Ising The problem of finding an approximate solution that is as small as possible (i.e., an approximate solution whose value of the objective function is as close as possible to the optimal value) is called the Ising problem.

[0061] Since the spin s of formula (1) i is a binary variable, so by using the formula (1+s i ) / 2, it is easy to transform the discrete variables (bits) used in the combinatorial optimization problem. Therefore, by transforming the combinatorial optimization problem into the Ising problem and making the Ising machine perform the calculation, the solution to the combinatorial optimization problem can be found. The problem of finding a solution that minimizes a quadratic objective function with discrete variables (bits) that take a certain value of 0 or 1 as variables is called the QUBO (Quadratic Unconstrained Binary Optimization) problem. The Ising problem expressed in formula (1) can be said to be equivalent to the QUBO problem.

[0062] For example, quantum annealing machines, coherent Ising machines, or quantum branching machines have been proposed as hardware installations of the Ising machine. The quantum annealing machine uses a superconducting circuit to implement quantum annealing. The coherent Ising machine uses the oscillation phenomenon of a network formed by an optical parametric oscillator. The quantum branching machine uses the branching phenomenon of quantum mechanics in a network of parametric oscillators with the Kerr effect. These hardware installations have the potential to significantly shorten the computing time, but on the other hand, there are technical problems that make it difficult to achieve large-scale and stable operation.

[0063] Therefore, it is also possible to use widely used digital computers to solve the Ising problem. Compared with the hardware installation using the above-mentioned physical phenomenon, digital computers are easy to achieve large-scale and stable operation. As an example of an algorithm used to solve the Ising problem with a digital computer, simulated annealing (SA) can be cited. Technology for performing simulated annealing at a higher speed is being developed. However, since conventional simulated annealing is a sequential update algorithm in which each variable is updated sequentially, it is difficult to achieve high-speed calculation processing due to parallelization.

[0064] In view of the above technical problems, a simulation branching algorithm is proposed for solving large-scale combinatorial optimization problems at high speed by parallel calculations of digital computers. Hereinafter, an information processing device, an information processing system, an information processing method, a storage medium and a program for solving combinatorial optimization problems using the simulation branching algorithm are described.

[0065] First, an overview of the simulation branching algorithm is described.

[0066] In the simulation branching algorithm, for two variables x with N i ,y i (i=1, 2, ..., N), numerically solve the following simultaneous ordinary differential equations (2). N variables x i Corresponding to the spin s of the Ising model i On the other hand, N variables y i They are equivalent to the amount of exercise respectively. Let the variable x i ,y i are all continuous variables. i A vector with (i=1, 2, ..., N) as its elements is called the first vector. i A vector whose elements are (i=1, 2, ..., N) is called a second vector.

[0067] [Formula 2]

[0068]

[0069] Here, H is the Hamiltonian of the following equation (3).

[0070] [Formula 3]

[0071]

[0072] In (2), instead of the Hamiltonian H in equation (3), the potential (constraint potential function) term G(x1, x2, ... x N The function including not only Hamiltonian H but also constraint potential function is called extended Hamiltonian, which is distinguished from the original Hamiltonian H.

[0073] [Formula 4]

[0074]

[0075] For example, it is possible to use a function G(x1, x2, ... x N ). However, the function G(x1, x2, ...xN ) is derived. In addition, in equation (4), the original Hamiltonian H is added with the function G (x1, x2, ... x N ). However, the extended Hamiltonian is not limited to include the function G(x1, x2, ...x N ) method.

[0076] By calculating the time advancement of the simulated branching algorithm, the variable x can be i ,y i The value of (i=1, 2, ..., N) is repeatedly updated to find the spin s of the Ising model i (i=1, 2, ..., N). That is, when performing time-advancing calculations, the values ​​of the elements of the first vector and the values ​​of the elements of the second vector are repeatedly updated. Here, each coefficient is explained assuming that time-advancing calculations are performed. However, the calculations of the analog branching algorithm can also be performed by methods other than time-advancing.

[0077] In (2) and (3), coefficient D corresponds to demodulation. Coefficient p(t) corresponds to pumping amplitude. In the case of time advancement of the calculation analog branching algorithm, coefficient p(t) increases monotonically with the number of updates. The initial value of coefficient p(t) can also be set to 0. Coefficient K corresponds to a positive Kerr coefficient. As coefficient c, a constant coefficient can be used. For example, the value of coefficient c can also be determined before the execution of the calculation based on the analog branching algorithm. For example, coefficient c can be set to be close to J. (2) The value of the reciprocal of the maximum eigenvalue of the matrix. For example, the value of c = 0.5D√(N / 2n) can be used. Here, n is the number of edges of the curve related to the combinatorial optimization problem. In addition, a(t) is a coefficient that increases with p(t) in the calculation of time advancement. For example, √(p(t) / K) can be used as a(t). In addition, the vector h of the local magnetic field in (3) and (4) is i Can be omitted.

[0078] When the value of the coefficient p(t) exceeds a predetermined value in the time-advanced calculation, if the variable x is a positive value in the first vector, i Transformed to +1, the variable x is negative i Transformed to -1, we can get i is a solution vector of the elements. This solution vector is equivalent to the solution of the Ising problem. In addition, the above-mentioned transformation can also be performed based on the number of updates of the first vector and the second vector to determine whether the solution vector is obtained.

[0079] When calculating the time progression of the analog branching algorithm, the Symplectic Euler method can be used to transform the above (2) into a discrete recursive equation and solve it. The following (5) shows an example of the analog branching algorithm after the recursive equation is transformed.

[0080] [Formula 5]

[0081] x i (t+Δt)=x i (t)+Dy i (t)Δt (5)

[0082]

[0083] Here, t is the time and Δt is the time step (time scale width). In addition, in (5), the time t and the time step Δt are used to express the corresponding relationship with the differential equation. However, when the algorithm is actually installed in software or hardware, it is not necessary to include the time t and the time step Δt as explicit parameters. For example, if the time step Δt is set to 1, the time step Δt can be removed from the algorithm during installation. In the case where the time t is not included as an explicit parameter when installing the algorithm, just replace x in (4). i (t+Δt) is interpreted as x i That is, "t" in the above (4) represents the value of the variable before the update, and "t+Δt" represents the value of the variable after the update.

[0084] When calculating the time progress of the simulated branching algorithm, based on the variable x after increasing the value of p(t) from an initial value (for example, 0) to a predetermined value, i The sign of the spin s can be found i For example, if you use x i >0 is sgn(x i )=+1, at x i <0 is sgn(x i )=-1, then when the value of p(t) increases to a specified value, by changing the variable x i By using the symbolic function transformation, we can find the spin s i As a symbolic function, for example, we can use i ≠0 is sgn(x i )=x i / |x i |、When x i = 0 is sgn(x i )=+1 or -1. Regarding the solution of the combinatorial optimization problem (for example, the spin s of the Ising modeli ) is not particularly limited. For example, the solution (solution vector) of the combinatorial optimization problem may be obtained when the number of updates of the first vector and the second vector or the value of the first coefficient p becomes larger than a threshold value.

[0085] In this way, the processing circuit of the information processing device may be configured to calculate the solution by transforming the first variable, which is a positive value, into the first value and transforming the first variable, which is a negative value, into a second value that is smaller than the first value. In addition, the processing circuit may be configured to determine whether to calculate the solution based on the value of the first coefficient or the number of updates of the first vector and the second vector. Here, the first value is, for example, +1. The second value is, for example, -1. However, the first value and the second value may also be other values.

[0086] Figure 6 The flowchart of shows an example of processing when the solution of the simulation branching algorithm is calculated by time advancement. Figure 6 While explaining the processing.

[0087] First, the computing server obtains the matrix J corresponding to the problem from the management server 1 ij and vector h i (Step S101). Next, the computing server initializes the coefficients p(t) and a(t) (Step S102). For example, in step S102, the coefficients p and a can be set to 0, but the initial values ​​of the coefficients p and a are not limited. Next, the computing server initializes the first variable x i and the second variable y i Initialization (step S103). Here, the first variable x i is the element of the first vector. In addition, the second variable y i is an element of the second vector. In step S103, the computing server may also use a pseudo-random number to replace x i and i Initialization. However, it does not limit x i and i In addition, variables may be initialized at a different timing, or at least a certain variable may be initialized multiple times.

[0088] Next, the calculation server calculates the corresponding second vector element y by i Weighted and combined with the first vector element x i Add and update the first vector (step S104). For example, in step S104, the variable x i Add Δt×D×y i Furthermore, the computing server converts the element y of the second vector iUpdate (steps S105 and S106). For example, in step S105, the variable y i Add Δt×[(p-D-K×x i × i )×x i In step S106, the variable y can be i Add -Δt×c×h i ×a-Δt×c×ΣJ ij × j .

[0089] Next, the computing server updates the values ​​of coefficients p and a (step S107). For example, a certain value (Δp) is added to coefficient p, and coefficient a is set to the positive square root of the updated coefficient p. However, as described later, this is only one example of a method for updating the values ​​of coefficients p and a. Next, the computing server determines whether the number of updates of the first vector and the second vector is less than a threshold value (step S108). When the number of updates is less than the threshold value (“Yes” in step S108), the computing server executes the processing of steps S104 to S107 again. When the number of updates is greater than the threshold value (“No” in step S108), based on the element x of the first vector, the calculation server performs the processing of steps S104 to S107 again. i Find the spin s of the elements of the solution vector i (Step S109). In step S109, for example, in the first vector, the variables x and i Transformed to +1, the variable x is negative i Transform it to -1 and get the solution vector.

[0090] In addition, in the determination of step S108, when the number of updates is less than the threshold value ("Yes" in step S108), the value of the Hamiltonian may be calculated based on the first vector, and the first vector and the value of the Hamiltonian may be stored. Thus, the user can select the approximate solution closest to the optimal solution from a plurality of first vectors.

[0091] Alternatively, you can Figure 6 At least one of the processes shown in the flowchart of is executed in parallel. For example, the processes of steps S104 to S106 may be executed in parallel to update at least a portion of the N elements respectively possessed by the first vector and the second vector in parallel. For example, the processes may be parallelized using multiple computing servers. The processes may also be parallelized using multiple processors. However, the installation used to realize the parallelization of the processes and the method of the parallelization of the processes are not limited.

[0092] The variable x represented in the above steps S105 to S106 is i and iThe execution order of the update processing of is just an example. Therefore, the variable x can also be executed in a different order. i and i For example, you can also execute the variable x i Update process and variable y i In addition, the order of sub-processes included in the update process of each variable is not limited. For example, variable y i The execution order of the addition processing included in the update processing can also be the same as Figure 6 The execution order and timing of the processing that is the prerequisite for executing the updating processing of each variable are not particularly limited. For example, the calculation processing of the problem item can be combined with the processing including the variable x. i The update process of the variable x is executed in parallel with the other processes. i and i The order and timing in which the update processing, the sub-processing included in the update processing of each variable, and the calculation processing of the question item are executed are not limited, and the same applies to the processing of each flowchart shown later.

[0093] [Objective function including constraints]

[0094] As described above, instead of the Hamiltonian H, a constraint potential function G(x1, x2, ... x N ) is used to simulate the branching algorithm. N ) is described in detail.

[0095] Function G(x1, x2, ...x N ) can also be a function representing a constraint. In addition, as described below, the function G(x1, x2, ... x N ) can also be a function representing multiple constraints. Here, let g m (x1, x2, ... x N )(m=1, 2...M) represents a function (constraint function) that expresses M constraints respectively. In this case, multiple functions g can be used. m (x1, x2, ... x N ) defines the above function G(x1, x2, ...x N ).

[0096] In the calculation of the optimization problem including the constraint condition, it is preferable to use the function G (x1, x2, ... x N ). Therefore, the constraint condition is preferably as follows: m (x1, x2, ... x N)=0(m=1,2…M). However, it is not limited to expressing all constraint conditions as equality constraints. Therefore, the constraint conditions may be transformed into constraint functions based on the following rules (a) to (c). Here, g * m (x1, x2, ... x N ) is the constraint function before transformation.

[0097] (a) In such m * (x1, x2, ... x N )=0, if the constraint condition is originally an equality constraint, then g m (x1, x2, ... x N ) = g m * (x1, x2, ... x N ). Figure 7 The constraint x1+x2=0 represents an example of an equality constraint. Figure 7 The table represents the condition table under the constraint x1+x2=0. Figure 7 The rows in the table that are not grayed out correspond to the combinations of variables that satisfy the constraints.

[0098] (b) When the original constraint is g m * (x1, x2, ... x N )≤0, we can set g m (x1, x2, ... x N )=max(0,g m * (x1, x2, ... x N )). Here, the function max is a function that takes the larger of the first and second independent variables as its value. Figure 8 The constraint x1+x2≤0 is an example of an inequality constraint. Figure 8 The table represents the condition table of the constraint x1+x2≤0. Figure 8 The rows in the table that are not grayed out correspond to the combinations of variables that satisfy the constraints.

[0099] (c) Under the original constraint condition g m * (x1, x2, ... x N )≥0, we can set g m (x1, x2, ... x N )=min(0,g m * (x1, x2, ... x N)). Here, the function min is a function that takes the smaller of the first and second arguments as its value. Fig. 9 The constraint x1+x2≥0 represents an example of an inequality constraint. Fig. 9 The table represents the condition table of the constraint x1+x2≥0. Fig. 9 The non-greyed-out rows in the table correspond to the combinations of variables that satisfy the constraints.

[0100] If the above-mentioned transformation rules (a) to (c) are applied to the constraint conditions, the extended Hamiltonian H′ represented by the following equation (6) can be used.

[0101] [Formula 6]

[0102]

[0103] c in formula (6) g (t) is, for example, a coefficient that increases monotonically according to the number of updates.

[0104] In this case, the effect of the value of the function (7) described below becomes larger as the number of updates increases.

[0105] [Formula 7]

[0106]

[0107] When the objective function of equation (6) is used, the simultaneous ordinary differential equations represented by the following equation (8) are executed for each of the two variables x having N numbers. i ,y i (i=1, 2, ..., N) is a process of numerically solving.

[0108] [Formula 8]

[0109]

[0110] Using the symplectic Euler method, the above (8) can be transformed into a discrete recursive formula and the calculation of the analog branching algorithm can be performed. The following (9) shows an example of the analog branching algorithm after the recursive formula is transformed.

[0111] [Formula 9]

[0112] x i (t+Δt)=x i (t)+Dy i (t)Δt (9)

[0113]

[0114] The term of (10) below in (9) is derived from the Ising energy. The form of this term is determined by the problem to be solved, so it is called a problem term.

[0115] [Formula 10]

[0116]

[0117] On the other hand, the following term (11) in (9) corresponds to the variable y i The first amendment to

[0118] [Formula 11]

[0119]

[0120] As the coefficient c g The initial value of (t) can be a relatively small value (0 or near 0). i The initial value set in (i = 1, 2, ..., N) can prevent the value of the term in (7) from becoming too large, which may cause the stability of the calculation process to be impaired. In addition, not only the vicinity of a specific local solution, but also a larger domain can be explored, thereby increasing the probability of finding an optimal solution or an approximate solution close to it. Furthermore, since it is no longer necessary to set the time step Δt of the algorithm in (9) to a small value, the calculation time can be suppressed.

[0121] For example, the coefficient c defined based on the following equation (12) can be used: g (t).

[0122] [Formula 12]

[0123] c g (t) = c g (0)+Δc g t (12)

[0124] Here, c g (0) is the coefficient c g (t) initial value. Δc g is the coefficient by which the update count or counting variable t is multiplied. g (0) and Δc g Positive values ​​can be used. If the coefficient c defined in (12) is used g (t), the value of the correction term in (11) increases according to the number of updates.

[0125] The coefficient c explained here g The initial value of (t) and the coefficient c corresponding to the number of updates g The variation of (t) is just an example. Therefore, the coefficient cg The initial value of (t) and the coefficient c corresponding to the number of updates g The variation pattern of (t) may also be different from the above.

[0126] The processing circuit of the information processing device may be configured to calculate a function (constraint term) including the product of the constraint function and the derivative obtained by differentiating the constraint function with respect to a certain element (first variable) of the first vector. In addition, the processing circuit may be configured to calculate the above-mentioned product for each of a plurality of constraint conditions and to calculate the constraint term by adding the plurality of products. Here, the above-mentioned function g m is an example of a constraint function.

[0127] Fig.10 The flowchart of FIG. 1 shows an example of the solution process of the extended Hamiltonian including the first correction term. Fig.10 While explaining the processing.

[0128] First, the computing server obtains the matrix J corresponding to the problem from the management server 1 ij and vector h i (Step S111). The computing server then initializes the coefficients p(t) and a(t) (Step S112). For example, in step S112, the coefficients p and a may be set to 0, but the initial values ​​of the coefficients p and a are not limited. Next, the computing server calculates the value of the coefficients p(t) and a(t) based on the corresponding element y of the second vector. i The first element x of the vector i Update (step S113). For example, in step S113, the variable x i Add Δt×D×y i Furthermore, the computing server converts the element y of the second vector i Update (steps S114 to S116). For example, in step S114, the variable y i Add Δt×[(p-D-K×x i × i )×x i In step S115, the variable y can be i Add -Δt×c×h i ×a-Δt×c×ΣJ ij × j In step S116, the variable y can be i Plus

[0129] Next, the computing server calculates the coefficients p and c g For example, a certain value (Δp) can be added to the coefficient p, the coefficient a can be set to the positive square root of the updated coefficient p, and the coefficient c can be set to the positive square root of the updated coefficient p.g Add Δc g t. However, this is just the coefficients p, c g and a. In addition, the computing server determines whether the number of updates of the first vector and the second vector is less than a threshold value (step S118). If the number of updates is less than the threshold value (“Yes” in step S118), the computing server repeatedly executes the processing of steps S113 to S117. If the number of updates is greater than the threshold value (“No” in step S118), the computing server repeatedly executes the processing of steps S113 to S117. i Find the spin s of the elements of the solution vector i (Step S119) In step S119, for example, in the first vector, the variables x and i Transformed to +1, the variable x that is negative i Transform it to -1 and get the solution vector.

[0130] In addition, the computing server can also calculate the value of the objective function based on the first vector and the second vector at any timing. The computing server can save the first vector, the second vector and the value of the objective function in the storage unit. These processes can also be performed each time when the judgment of step S118 is affirmative. In addition, these processes can also be performed at a timing that is part of the timing when the judgment of step S118 becomes affirmative. Furthermore, the above-mentioned processes can also be performed at other timings. The user can determine the frequency of calculating the value of the objective function based on the available storage area and the amount of computing resources. It is also possible to determine whether to continue the loop processing based on whether the number of combinations of the first vector, the second vector and the value of the objective function stored in the storage unit exceeds a threshold at the timing of step S118. In this way, the user can select the first vector closest to the optimal solution from multiple first vectors stored in the storage unit to calculate the solution vector.

[0131] Alternatively, you can Fig.10 At least one of the processes shown in the flowchart may be executed in parallel. For example, the processes of steps S113 to S116 may be executed in parallel to update at least a portion of the N elements respectively possessed by the first vector and the second vector in parallel. For example, the processes may be parallelized using multiple computing servers. The processes may also be parallelized by multiple processors. However, the installation and method of the parallelization of the processes are not limited.

[0132] By executing Fig.10 The processing shown in the flowchart can obtain the optimal solution that satisfies the constraints or an approximate solution close to it in a relatively short time.

[0133] [Calculation based on the extended Lagrangian method]

[0134] The above-mentioned formula (6) is only an example of an extended Hamiltonian that can reflect the constraint condition. For example, as in the following formula (13), further terms may be added to the extended Hamiltonian for stable calculation.

[0135] [Formula 13]

[0136]

[0137] The extended Hamiltonian H″ shown in equation (13) includes both the penalty function (the first term of the second stage) and the Lagrangian function (the second term of the second stage). This method of including both the penalty function and the Lagrangian function in the extended Hamiltonian is called the extended Lagrangian method.

[0138] When the extended Hamiltonian of equation (13) is used, the simultaneous ordinary differential equations represented by the following equation (14) are executed for two variables x each having N number of variables. i ,y i (i=1, 2, ..., N) is a process of solving the problem numerically.

[0139] [Formula 14]

[0140]

[0141] In addition, the coefficient λ can be used in the calculation of the simultaneous ordinary differential equations (14). m Update based on the following formula (15).

[0142] [Formula 15]

[0143]

[0144] The vector x in formula (15) min is a vector equivalent to a local solution (x1, x2, ..., x N ). The local solution x min For example, it can be obtained by applying a local heuristic or best-first heuristic algorithm to the first vector in the calculation. As an example of a local heuristic, a negative hill climbing algorithm can be cited. When a negative hill climbing algorithm is used, for example, for each variable x i (i=1, 2, ..., N) Partially differentiate the extended Hamiltonian H″ and obtain the following partial derivatives of (16).

[0145] [Formula 16]

[0146]

[0147] Furthermore, based on each partial derivative, a vector whose evaluation value is smaller than that of the extended Hamiltonian H″ is searched in the vicinity of the first vector. For example, a vector obtained by adding the product of the partial derivative of (16) and Δx to each element of the first vector can be used in the next iteration. In the next iteration, the above-mentioned partial derivatives are also calculated, and a vector whose evaluation value of the extended Hamiltonian H″ becomes smaller is searched in the vicinity of the vector. This process is repeated until the values ​​of the partial derivatives of (16) are approximately zero. For example, the absolute value of the partial derivative of (16) can be compared with a threshold value to determine whether to continue the iterative process. The iteratively processed vector (x1, x2, ..., x N ) as a local solution. The local solution exploration method described here is just an example. Therefore, other algorithms can also be used to explore local solutions.

[0148] Using the symplectic Euler method, the above (14) and (15) can be transformed into discrete recursive expressions, and the calculation of the analog branching algorithm can be performed. The following (17) shows an example of the analog branching algorithm after the transformation into the recursive expression.

[0149] [Formula 17]

[0150] x i (t+Δt)=x i (t)+Dy i (t)Δt (17)

[0151]

[0152] λ m (t+Δt)=λ m (t)+Δtc g (t)g m (x min )

[0153] In the algorithm (17), the coefficient c g For example, the coefficient c can be updated based on the above formula (12). g renew.

[0154] When the constraints are not met, g m Therefore, during the period when the constraint condition is not satisfied, the coefficient λ of the Lagrangian term in (17) m The rate of increase becomes larger. Therefore, the effect of the Lagrangian term becomes larger, and the first vector and the second vector can be changed in the direction that satisfies the constraint condition. Therefore, if the algorithm (17) including the Lagrangian function is used, the coefficient c can be g The initial value of is set smaller or the coefficient c gThe value of increases more gradually. This can prevent the slope of the potential of the extended Hamiltonian from becoming too large, and can stabilize the calculation process.

[0155] The processing circuit of the information processing device may be configured to calculate a second correction term including a product of a third coefficient and a derivative obtained by differentiating the constraint condition with respect to a first variable, and to add the second correction term to the second variable. In addition, the processing circuit may be configured to calculate the above-mentioned product for each of a plurality of constraint conditions, and to calculate the second correction term by adding the plurality of products. Here, the above-mentioned coefficient λ m It is an example of the third coefficient. The above-mentioned Lagrangian term is an example of the second correction term.

[0156] In addition, the processing circuit may be configured to increase the absolute value of the third coefficient according to the number of updates. Furthermore, the processing circuit may be configured to calculate the evaluation value of the constraint condition by substituting the elements of the first vector corresponding to the local solution of the objective function (extended Hamiltonian) into the constraint condition, and add the third coefficient to the product of the second coefficient and the evaluation value. For example, the vector x calculated by the local heuristic method can be min Substitute into the function g m , and obtain the evaluation value.

[0157] Fig.11 The flowchart of shows an example of the solution process of the extended Hamiltonian which also includes the second correction term. Fig.11 While explaining the processing.

[0158] First, the computing server initializes the coefficients p(t) and a(t) (step S121). For example, in step S121, the values ​​of the coefficients p and a can be set to 0, but the initial values ​​of the coefficients p and a are not limited. Fig.11 It is not shown in the figure, but it is assumed that the calculation server obtains the matrix J corresponding to the problem from the management server 1. ij and vector h i Next, the server calculates the element y based on the corresponding second vector i The first element x of the vector i Update (step S122). For example, in step S122, the variable x i Add Δt×D×y i .

[0159] Then, the computing server converts the element y of the second vector i Update (steps S123 to S126). For example, in step S123, the variable y i Add Δt×[(p-D-K×x i × i )×x iIn step S124, the variable y can be i Add -Δt×c×h i ×a-Δt×c×ΣJ ij × j In step S125, the variable y can be i Plus In step S126, the variable y may be i Plus

[0160] Then, the computing server calculates the local solution of the extended Hamiltonian H″ (step S127). In step S127, the local solution can be calculated by the negative hill climbing algorithm as described above. However, the local solution can also be calculated by other algorithms. Next, the coefficient λ is calculated based on the local solution calculated in the previous step. m Update (step S128). For example, in step S128, the local solution can be substituted into the function g m In the equation, find the function g m After the value of m Add Δt×c g ×g m .

[0161] Next, the computing server calculates the coefficients p and c g For example, a certain value (Δp) may be added to the coefficient p, the coefficient a may be set to the positive square root of the updated coefficient p, and the coefficient c may be set to the positive square root of the updated coefficient p. g Add Δc g t. However, this is just the coefficients p, c g and a. In addition, the computing server determines whether the number of updates of the first vector and the second vector is less than a threshold value (step S130). If the number of updates is less than the threshold value (“Yes” in step S130), the computing server repeatedly executes the processing of steps S122 to S129. If the number of updates is greater than the threshold value (“No” in step S130), the computing server repeatedly executes the processing of steps S122 to S129. i Find the spin s of the elements of the solution vector i (Step S131). In step S131, for example, in the first vector, the variables x and i Transformed to +1, the variable x is negative i Transform it to -1 and get the solution vector.

[0162] In addition, the computing server may also calculate the value of the objective function based on the first vector and the second vector at any timing. The computing server may save the first vector, the second vector and the value of the objective function in the storage unit. These processes may also be performed each time when the determination in step S130 is affirmative. In addition, these processes may also be performed at a portion of the timings when the determination in step S130 is affirmative. Furthermore, the above-mentioned processes may also be performed at other timings. The user can calculate the frequency of the value of the objective function based on the available storage area and the amount of computing resources. It is also possible to determine whether to continue the loop processing based on whether the number of combinations of the first vector, the second vector and the value of the objective function stored in the storage unit exceeds a threshold value at the timing of step S130. In this way, the user can select the first vector closest to the optimal solution from the multiple first vectors stored in the storage unit to calculate the solution vector.

[0163] Alternatively, you can Fig.11 At least one of the processes shown in the flowchart may be executed in parallel. For example, the processes of steps S122 to S126 may be executed in parallel to update at least a portion of the N elements of the first vector and the second vector in parallel. For example, the processes may be parallelized using multiple computing servers. The processes may also be parallelized by multiple processors. However, the installation and method of the parallelization of the processes are not limited.

[0164] By executing Fig.11 The processing shown in the flowchart can obtain the optimal solution that satisfies the constraints or an approximate solution close to it in a relatively short time.

[0165] [Reduction in calculation amount]

[0166] In the above Fig.11 In the flowchart shown, each time the coefficient λ m However, in the information processing device of this embodiment, the coefficient λ may not necessarily be updated every time. m Update. For example, the coefficient λ can be updated in some iterations. m Updates skipped.

[0167] For example, instead of the above (15), the coefficient λ may be changed based on the following (18): m For example, the coefficient λ may be updated based on the rule (18) below for each W times (W is a positive integer). m Updates.

[0168] [Formula 18]

[0169] λ m (t+Δt)=λ m (t)+ΔtZc g(t)g m (x min ) (18)

[0170] The z in (18) can be a fixed value or a variable. For example, a value equal to W can be set as Z. In addition, 1 / (c g ×g m ).

[0171] Execution coefficient λ m The frequency of the update process can also be changed. For example, the value of W can be changed to 1 / (c g ×g m ). Thus, the constraint condition is satisfied, and the smaller the value of g, the smaller the coefficient λ m The smaller the update frequency, the more computational effort can be suppressed.

[0172] In the coefficient λ m In the iteration in which the update process of is skipped, the calculation process of the local solution of the extended Hamiltonian H″ can also be skipped. In this way, the required amount of calculation can be reduced while maintaining the stability of the above-mentioned calculation process.

[0173] In this way, the processing circuit of the information processing device may be configured to increase the absolute value of the third coefficient in a part of the update times.

[0174] Fig.12 The flow chart represents the coefficient λ m An example of a solution process when part of the update process of is skipped. Fig.12 While explaining the processing.

[0175] First, the calculation server initializes the coefficients p(t) and a(t) (step S141). For example, the coefficients p and a can be set to 0 in step S141, but the initial values ​​of the coefficients p and a are not limited. The counting variable cnt described later can be initialized to 0. In addition, although Fig.12 It is not shown in the figure, but it is assumed that the calculation server obtains the matrix J corresponding to the problem from the management server 1. ij and vector h i Next, the server calculates the element y based on the corresponding second vector i The first element x of the vector i Update (step S142). For example, in step S142, the variable x i Add Δt×D×y i .

[0176] Then, the computing server converts the element y of the second vector iUpdate (steps S142 to S146). For example, in step S143, the variable y i Add Δt×[(p-D-K×x i × i )×x i In step S144, the variable y can be i Add -Δt×c×h i ×a-Δt×c×ΣJ ij × j In step S145, the variable y can be i Plus In step S146, the variable y may be i Plus

[0177] Then, the calculation server increments the count variable cnt and determines whether the count variable cnt is a multiple of W (S147). If the count variable cnt is not a multiple of W ("No" in S147), the calculation server increments the coefficients p and c g For example, a certain value (Δp) can be added to the coefficient p, the coefficient a can be set to the positive square root of the updated coefficient p, and the coefficient c can be set to the positive square root of the updated coefficient p. g Add Δc g t. However, this is just the coefficients p and c g And an example of a method for updating the value of a.

[0178] When the count variable cnt is a multiple of W (“Yes” in S147), the computing server calculates a local solution of the extended Hamiltonian H″ (step S148). In step S148, the local solution can be calculated by the negative hill climbing algorithm as described above. However, the local solution can also be calculated by other algorithms. Next, the coefficient λ is calculated based on the local solution calculated in the previous step. m Update (step S149). For example, in step S149, the local solution can be substituted into the function g m In finding the function g m After the value of coefficient λ m Add Δt×c g ×g m After the processing of step S149, the process proceeds to the processing of the above-mentioned step S150.

[0179] After the processing of step S150, the computing server determines whether the update times of the first vector and the second vector are less than a threshold value (step S151). If the update times are less than the threshold value (“Yes” in step S151), the computing server executes the processing after step S142 again. If the update times are greater than the threshold value (“No” in step S151), the computing server performs the processing after step S142 again. i Find the spin s of the elements of the solution vector i (Step S152) In step S152, for example, in the first vector, variables x and i Transformed to +1, the variable x that is negative i Transform it to -1 and get the solution vector.

[0180] In addition, the calculation server may also calculate the value of the objective function based on the first vector and the second vector at any timing. The calculation server may save the first vector, the second vector and the value of the objective function in the storage unit. These processes may also be performed each time when the determination in step S151 is affirmative. In addition, these processes may also be performed at a portion of the timings when the determination in step S151 is affirmative. Furthermore, the above-mentioned processes may also be performed at other timings. The user can determine the frequency of calculating the value of the objective function based on the available storage area and the amount of computing resources. It is also possible to determine whether to continue the loop processing based on whether the number of combinations of the first vector, the second vector and the value of the objective function stored in the storage unit exceeds a threshold value at the timing of step S151. In this way, the user can select the first vector closest to the optimal solution from the multiple first vectors stored in the storage unit to calculate the solution vector.

[0181] Alternatively, you can Fig.12 At least one of the processes shown in the flowchart of is executed in parallel. For example, the processes of steps S142 to S146 may be executed in parallel to update at least a portion of the N elements respectively possessed by the first vector and the second vector in parallel. For example, the processes may be parallelized using multiple computing servers. The processes may also be parallelized using multiple processors. However, the installation and method of the parallelization of the processes are not limited.

[0182] By executing Fig.12 The processing shown in the flowchart can obtain the optimal solution or an approximate solution close to the optimal solution that satisfies the constraints in a relatively short time while suppressing the amount of calculation.

[0183] Hereinafter, examples of the information processing system, information processing method, storage medium, and program according to the present embodiment will be described.

[0184] The information processing system may also include a storage device and an information processing device. The storage device is configured, for example, to store a first variable as an element of a first vector and a second variable as an element of a second vector. The information processing device is configured, for example, to update the first variable based on the corresponding second variable, weight the first variable with a first coefficient and add the corresponding second variable, calculate the problem item using a plurality of first variables, add the second variable to the problem item, calculate the first correction item including the product of the constraint item and the second coefficient, add the second variable to the first correction item, and increase the absolute values ​​of the first coefficient and the second coefficient according to the number of updates. Here, the constraint item may also be based on a constraint function representing a constraint condition and have the first variable as an independent variable.

[0185] Furthermore, the information processing system may include a plurality of information processing devices. Each of the information processing devices may be configured to update at least a portion of the first vector and at least a portion of the second vector in parallel.

[0186] In the information processing method, for example, a first vector having a first variable as an element and a second vector having a second variable corresponding to the first variable as an element are repeatedly updated. For example, the information processing method may also include: a step of updating the first variable based on the corresponding second variable; a step of weighting the first variable with a first coefficient and adding the corresponding second variable; a step of calculating a problem item using a plurality of first variables; a step of adding the problem item to the second variable; a step of calculating a constraint item having the first variable as an independent variable based on a constraint condition; a step of calculating a first correction item including the product of the second coefficient and the constraint item; a step of adding the first correction item to the second variable; and a step of increasing the absolute values ​​of the first coefficient and the second coefficient according to the number of updates. In addition, the program may also cause a computer to execute the steps of the above-mentioned information processing method. Furthermore, the storage medium may also be a non-temporary computer-readable storage medium storing the program.

[0187] [Calculation of terms including multi-body interactions]

[0188] By using the simulated branching algorithm, it is possible to solve combinatorial optimization problems with objective functions of degree 3 or higher. The problem of solving the combination of variables that minimizes the objective function of degree 3 or higher with binary variables is called the HOBO (Higher Order Binary Optimization) problem. When dealing with the HOBO problem, the following formula (19) can be used as the energy formula in the Ising model extended to higher orders.

[0189] [Formula 19]

[0190]

[0191] Here, J (n)is an n-th order tensor, which is the local magnetic field h that makes equation (1) i and the coupling coefficient matrix J. For example, the tensor J (1) Equivalent to the local magnetic field h i A vector of n-order tensor J (n) In the equation (15), when there are the same values ​​in multiple suffixes, the value of the element is 0. In equation (15), the terms up to the third order are expressed, but the terms higher than this can also be defined in the same way as equation (19). Equation (19) is equivalent to the energy of the Ising model including multi-body interaction.

[0192] In addition, both QUBO and HOBO can be said to be a type of unconstrained polynomial binary optimization (PUBO). That is, the combinatorial optimization problem with a quadratic objective function in PUBO is QUBO. In addition, the combinatorial optimization problem with a tertiary or higher objective function in PUBO can be said to be HOBO.

[0193] When using the analog branching algorithm to solve the HOBO problem, it is sufficient to replace the Hamiltonian H in the above equation (3) with the Hamiltonian H in the following equation (20).

[0194] [Formula 20]

[0195]

[0196] Furthermore, based on equation (20), the problem item represented by the following equation (21) is derived.

[0197] [Formula 21]

[0198]

[0199] (21) Problem item z i Take the second equation of (20) with respect to a variable x i (The first vector element) is partially differentiated. The variable x of the partial differentiation is i It varies according to the index i. Here, the variable x i The index i of is equivalent to the index specifying the element of the first vector and the element of the second vector.

[0200] When performing calculations including terms of multi-body interactions, the above recursive formula (17) is replaced with the following recursive formula (22).

[0201] [Formula 22]

[0202] x i (t+Δt)=x i(t)+Dy i (t)Δt (22)

[0203]

[0204] λ m (t+Δt)=λ m (t)+Δtc g (t)g m (x min )

[0205] (22) is equivalent to further generalizing the recursive formula of (17). Similarly, in the recursive formula of (9) above, the term of multi-body interaction can also be used.

[0206] The questions described above are merely examples of questions that can be used by the information processing device of this embodiment. Therefore, the form of the question used in the calculation may be different from these.

[0207] [Algorithm Variations]

[0208] Here, a modification example of the analog branching algorithm is described. For example, various modifications may be made to the above-mentioned analog branching algorithm for the purpose of reducing errors or shortening calculation time.

[0209] For example, in order to reduce the calculation error, additional processing may be performed when updating the elements of the first vector. i When the absolute value of becomes greater than 1 due to updating, the element x of the first vector is i The value of sgn(x i ). That is, after being updated to become x i >1, variable x i The value of is set to 1. In addition, when x is updated i <-1, variable x i The value of is set to -1. Thus, the variable x can be used i Will spin s i By including this processing, the algorithm is similar to i = ±1, the physical model of N particles with walls at positions becomes equivalent. More generally, the operation circuit may be configured to set the first variable having a value smaller than the second value to the second value and to set the first variable having a value larger than the first value to the first value.

[0210] Then, when x is updated i >1, it can also be used for variable x i The corresponding variable y iMultiply by the coefficient \(r_f\). For example, if the coefficient \(r_f\) with \(-1 \lt r \leq 0\) is used, the above-mentioned wall becomes a wall with a reflection coefficient \(r_f\). In particular, when the coefficient \(r_f = 0\) is used, the algorithm becomes equivalent to the physical model where there are walls that cause perfectly inelastic collisions at positions \(x\) i \(=\pm1\). More generally, the arithmetic circuit can also be configured to update the second variable corresponding to the first variable with a value smaller than the first value or the second variable corresponding to the first variable with a value larger than the second value to the value obtained by multiplying the original second variable by the second coefficient. For example, the arithmetic circuit can also be configured to update the second variable corresponding to the first variable with a value smaller than \(-1\) or the second variable corresponding to the first variable with a value larger than \(1\) to the value obtained by multiplying the original second variable by the second coefficient. Here, the second coefficient corresponds to the above-mentioned coefficient \(r_f\).

[0211] In addition, the arithmetic circuit can also, when \(x\) i \gt 1\) after the update, set the value of the variable \(y\) i corresponding to the variable \(x\) i to a pseudo-random number. For example, a random number in the range \([-0.1, 0.1]\) can be used. That is, the arithmetic circuit can also be configured to set the value of the second variable corresponding to the first variable with a value smaller than the second value or the value of the second variable corresponding to the first variable with a value larger than the first value to a pseudo-random number.

[0212] If the update process is performed as above to prevent \(|x\) i |\gt 1\), then even if the non-linear terms \(K\times x\) i 2 in (5), (9), and (22) are removed, the value of \(x\) i will not diverge. Therefore, the algorithm expressed in the following (23) can be used.

[0213] [Equation 23]

[0214] \(x\) i (t + \Delta t)=x i (t)+Dy i (t)\Delta t\ (23)

[0215]

[0216] \(\lambda\) m (t + \Delta t)=\lambda m (t)+\Delta tc g (t)g m (x min )

[0217] In the algorithm of (23), a continuous variable x is used in the problem term instead of a discrete variable. Therefore, an error may occur with the discrete variable used in the original combinatorial optimization problem. In order to reduce the error, as shown in the following (24), in the calculation of the problem term, the value sgn(x) obtained by transforming the continuous variable x with a sign function may be used instead of the continuous variable x.

[0218] [Formula 24]

[0219] x i (t+Δt)=x i (t)+Dy i (t)Δt (24)

[0220]

[0221] λ m (t+Δt)=λ m (t)+Δtc g (t)g m (x min )

[0222] In (24), sgn(x) is equivalent to spin s.

[0223] In (24), the coefficient α of the term containing the first-order tensor in the problem term can also be set to a constant (for example, α = 1). In the algorithm of (24), since the product of the spins appearing in the problem term must take a value of -1 or 1, when dealing with the HOMO problem with a high-order objective function, it is possible to prevent the occurrence of errors caused by product operations. As in the algorithm of (24) above, the data calculated by the computing server can also include the variable s i (i=1, 2, ..., N) is the spin vector of the element (s1, s2, ..., s N ). By transforming the elements of the first vector using a sign function, a spin vector can be obtained.

[0224] [Example of parallelizing variable update processing]

[0225] Hereinafter, an example of parallelization of the update process of variables during calculation of the simulation branching algorithm will be described.

[0226] First, an example of installing a simulation branching algorithm on a PC cluster is described. The PC cluster is a system that connects multiple computers to achieve computing performance that cannot be obtained with a single computer. For example, Figure 1The information processing system 100 shown in FIG. 1 includes a plurality of computing servers and processors and can be used as a PC cluster. For example, in a PC cluster, by using MPI (Message Passing Interface), even a configuration in which memories are distributed among a plurality of computing servers such as the information processing system 100 can perform parallel computing. For example, the control program 14E of the management server 1, the computing program 34B of each computing server, and the control program 34C can be installed using MPI.

[0227] When the number of processors used by the PC cluster is Q, each processor can be made to perform the first vector (x1, x2, ..., x N ) contains the variable x i Similarly, each processor can perform the calculation of the second vector (y1, y2, ..., y N ) contains the variable y i That is, processor #j (j=1, 2, ..., Q) performs the calculation of the variables {x m |m=(j-1)L+1, (j-1)L+2,…,jL} and {y m |m=(j-1)L+1,(j-1)L+2,…,jL}. In addition, it is assumed that {y m The tensor J represented by the following (25) is required for the calculation of |m=(j-1)L+1,(j-1)L+2,…,jL} (n) It is stored in a storage area (eg, register, cache, memory, etc.) accessible to processor #j.

[0228] [Formula 25]

[0229]

[0230] Here, the case where each processor calculates a constant number of variables of the first vector and the second vector is described. However, the number of variables of the first vector and the second vector to be calculated may be different depending on the processor. For example, when there is a performance difference depending on the processor installed in the calculation server, the number of variables to be calculated may be determined according to the performance of the processor.

[0231] In order to change the variable y i The value of , needs the first vector (x1, x2, ..., x N ) values ​​of all components. The conversion to binary variables can be performed by using the sign function sgn(), for example. Therefore, the Allgather function can be used to make Q processors share the first vector (x1, x2, ..., x N) of all the components. Although for the first vector (x1, x2, ..., x N ), the value needs to be shared between processors, but for the second vector (y1, y2, ..., y N ) and tensor J (n) , it is not necessary to share values ​​between processors. Sharing of data between processors can be achieved, for example, by using inter-processor communication or storing data in a shared memory.

[0232] Processor #j calculates the problem item {z m |m=(j-1)L+1,(j-1)L+2,…,jL}. Then, processor #j calculates the value of the problem item {{z m |m=(j-1)L+1,(j-1)L+2,…,jL}, and change the variable {y m |m=(j-1)L+1,(j-1)L+2,…,jL} is updated.

[0233] As shown in the above equations, in the vector of the problem item (z1, z2, ..., z N ), the calculation of the tensor J(n) and the vector (x1, x2, ..., x N ) is the most computationally intensive process in the above algorithm and may become a bottleneck in improving the calculation speed. Therefore, in the installation of a PC cluster, the sum-of-products operation can be distributed to Q = N / L processors and executed in parallel, thereby shortening the calculation time.

[0234] Fig.13 An example of a multi-processor configuration is schematically shown. Fig.13 For example, the multiple computing nodes correspond to multiple computing servers of the information processing system 100. Fig.13 The high-speed link corresponds to, for example, an intermediate connection between computing servers formed by the cables 4 a to 4 c and the switch 5 of the information processing system 100 . Fig.13 The shared memory corresponds to the shared memory 32, for example. Fig.13 The processors of are, for example, equivalent to the processors 33A to 33D of each computing server. Fig.13 Although multiple computing nodes are shown in FIG. 1 , this does not prevent the use of a single computing node configuration.

[0235] exist Fig.13 The data configured in each component and the data transferred between components are shown in FIG. In each processor, the variable x is calculated. i ,y i In addition, the variable x is transferred between the processor and the shared memory. iIn the shared memory of each computing node, for example, the first vector (x1, x2, ..., x N ), the second vector (y1, y2, ..., y N ) of L variables and tensor J (n) Furthermore, in the high-speed link connecting the computing nodes, for example, the first vector (x1, x2, ..., x N ). When using the Allgather function, each processor will i Update, need the first vector (x1, x2, ..., x N )’s full elements.

[0236] in addition, Fig.13 The arrangement and transfer of data shown in FIG. 1 are merely examples. There are no particular limitations on the method of arranging data, the method of transferring data, and the method of implementing parallelization in the PC cluster.

[0237] In addition, GPU (Graphics Processing Unit) can also be used to perform calculations of simulated branching algorithms.

[0238] Fig.14 An example of a configuration using a GPU is schematically shown. Fig.14 In , multiple GPUs are connected to each other by high-speed links. Each GPU is equipped with multiple cores that can access a shared memory. Fig.14 In the configuration example, multiple GPUs are connected via a high-speed link to form a GPU cluster. Figure 1 In the case of each computing server, the high-speed link corresponds to the intermediate connection between the computing servers formed by the cables 4a to 4c and the switch 5. Fig.14 In the configuration example, multiple GPUs are used, but parallel calculations can also be performed using a single GPU. That is, Fig.14 Each GPU can perform Fig.13 That is, the processor of the information processing device (computing server) can also be the core of the Graphics Processing Unit (GPU).

[0239] In the GPU, the variable x i and i And the tensor J (n) is defined as a device variable. The GPU can calculate the variables y in parallel through the matrix-vector product function. i The tensor J required in the update (n) with the first vector (x1, x2, ..., x N). In addition, by repeatedly performing the matrix-vector product operation, the product of the tensor and the vector can be obtained. In addition, for the first vector (x1, x2, ..., x N ) and the calculation of the second vector (y1, y2, ..., y N ) can enable each thread to execute the i-th element (x i ,y i ) to realize parallel processing.

[0240] [Universal processing for solving combinatorial optimization problems]

[0241] The following describes the overall processing performed to solve the combinatorial optimization problem using the simulated branching algorithm.

[0242] Fig.15 The flowchart of FIG. 1 shows an example of the overall processing performed to solve the combinatorial optimization problem. Fig.15 While explaining the processing.

[0243] First, the combinatorial optimization problem is formulated (step S201). Then, the formulated combinatorial optimization problem is transformed into an Ising problem (in the form of an Ising model) (step S202). Next, the solution to the Ising problem is calculated by an Ising machine (information processing device) (step S203). Then, the calculated solution is verified (step S204). For example, in step S204, it is confirmed whether the constraint conditions are satisfied. In addition, in step S204, it is also possible to refer to the value of the energy function (Hamiltonian) to confirm whether the obtained solution is the optimal solution or an approximate solution close to it.

[0244] Then, it is determined whether to perform recalculation based on at least one of the verification result or the number of calculations in step S204 (step S205). When it is determined that recalculation is to be performed ("Yes" in step S205), the processing of steps S203 and S204 is performed again. On the other hand, when it is determined that recalculation is not to be performed ("No" in step S205), a solution is selected (step S206). For example, in step S206, the selection can be made based on the satisfaction of the constraint conditions or at least one of the values ​​of the energy function. In addition, if multiple solutions are not calculated, the processing of step S206 can also be skipped. Finally, the selected solution is transformed into a solution to the combinatorial optimization problem, and the solution to the combinatorial optimization problem is output (step S207).

[0245] By using the information processing device, information processing system, information processing method, storage medium and computer program product described above, it is possible to calculate the solution of the combinatorial optimization problem in a practical time. As a result, the solution of the combinatorial optimization problem becomes easier, which can promote social innovation and scientific and technological progress.

[0246] In addition, the present invention is not limited to the above-mentioned embodiments as it is, and in the implementation stage, the constituent elements can be deformed and concretized within the scope of the main purpose. In addition, various inventions can be formed by appropriate combinations of multiple constituent elements disclosed in the above-mentioned embodiments. For example, several constituent elements can also be deleted from all constituent elements shown in the embodiments. Furthermore, constituent elements across different embodiments can also be appropriately combined.

[0247] Description of Reference Numerals

[0248] 1 Management Server

[0249] 2 Network

[0250] 3a, 3b, 3c Computing Servers

[0251] 4a, 4b, 4c cables

[0252] 5 Switch

[0253] 6 Client Terminal

[0254] 10 Processor

[0255] 11 Management Department

[0256] 12 Transformation

[0257] 13. Control Unit

[0258] 14 Storage

[0259] 14A Question Data

[0260] 14B Calculation Data

[0261] 14C Management Program

[0262] 14D Transformation Procedure

[0263] 14E, 34C Control Program

[0264] 15, 31 Communication circuit

[0265] 16 Input Circuit

[0266] 17 Output Circuit

[0267] 18 Operating device

[0268] 19 Display device

[0269] 20 Bus

[0270] 32 Shared Memory

[0271] 33A, 33B, 33C, 33D processors

[0272] 34 Storage Devices

[0273] 34A Calculation data

[0274] 34B Calculation Procedure

[0275] 35 Host Bus Adapter

Claims

1. An information processing device, wherein: have: a common storage unit configured to store a first variable and a second variable, the first variable being an element of a first vector and the second variable being an element of a second vector; and The processing circuit is configured to update the first variable based on the corresponding second variable, weight the first variable with a first coefficient and add the first variable to the corresponding second variable, calculate a problem item using a plurality of the first variables, add the second variable to the problem item, calculate a first correction item including a product of a constraint item and a second coefficient, add the second variable to the first correction item, and increase the absolute values ​​of the first coefficient and the second coefficient according to the number of updates, The constraint item is a constraint item based on a constraint function representing a constraint condition, and has the first variable as an independent variable. The above-mentioned processing circuit is multiple, Each of the processing circuits is configured to update at least a portion of the first vector and at least a portion of the second vector in parallel. The shared storage unit is shared by a plurality of the processing circuits. The first vector is transferred between the common storage unit and each of the plurality of processing circuits.

2. The information processing device according to claim 1, wherein: The processing circuit is configured to calculate the constraint term including a product of the constraint function and a derivative obtained by differentiating the constraint function with respect to the first variable.

3. The information processing device according to claim 2, wherein: The processing circuit is configured to calculate the product for each of the plurality of constraint conditions and to add the plurality of products to calculate the constraint term.

4. The information processing device according to claim 1, wherein: The processing circuit is configured to calculate a second correction term including a product of a third coefficient and a derivative obtained by differentiating the constraint condition with respect to the first variable, and to add the second correction term to the second variable.

5. The information processing device according to claim 4, wherein: The processing circuit is configured to calculate the product for each of the plurality of constraint conditions, and to add the plurality of products to calculate the second correction term.

6. The information processing device according to claim 4, wherein: The processing circuit is configured to increase the absolute value of the third coefficient according to the number of updates.

7. The information processing device according to claim 6, wherein: The processing circuit is configured to increase the absolute value of the third coefficient at a part of the update times.

8. The information processing device according to claim 6, wherein: The processing circuit is configured to calculate an evaluation value of the constraint function by substituting the first variable corresponding to a local solution of the objective function into the constraint function, and to add the third coefficient to a product of the second coefficient and the evaluation value.

9. The information processing device according to any one of claims 1 to 8, wherein: The above-mentioned problem items calculated by the above-mentioned processing circuit are based on the Ising model.

10. The information processing device according to claim 9, wherein: The above problem items calculated by the above processing circuit include multi-body interactions.

11. The information processing device according to any one of claims 1 to 8, wherein: The processing circuit is configured to calculate a solution by converting the first variable, which is a positive value, into a first value and converting the first variable, which is a negative value, into a second value smaller than the first value.

12. The information processing device according to claim 11, wherein: The processing circuit is configured to determine whether or not to calculate a solution based on the value of the first coefficient or the number of times the first vector and the second vector are updated.

13. An information processing system, wherein: have: a common storage unit configured to store a first variable and a second variable, the first variable being an element of a first vector and the second variable being an element of a second vector; and The information processing device is configured to update the first variable based on the corresponding second variable, weight the first variable with a first coefficient and add the first variable to the corresponding second variable, calculate a problem item using a plurality of the first variables, add the second variable to the problem item, calculate a first correction item including a product of a constraint item and a second coefficient, add the second variable to the first correction item, and increase the absolute values ​​of the first coefficient and the second coefficient according to the number of updates, The constraint item is a constraint item based on a constraint function representing a constraint condition, and has the first variable as an independent variable. The above-mentioned information processing device is multiple, Each of the information processing devices is configured to update at least a portion of the first vector and at least a portion of the second vector in parallel. The shared storage unit is shared by a plurality of the information processing devices. The first vector is transferred between the common storage unit and each of the plurality of information processing devices.

14. An information processing method, comprising: using a common storage unit and a plurality of processing circuits, wherein the common storage unit stores a first vector having a first variable as an element and a second vector having a second variable corresponding to the first variable as an element, and the plurality of processing circuits repeatedly updating the first vector having the first variable as an element and the second vector having the second variable corresponding to the first variable as an element, in, include: The step of each processing circuit updating the first variable based on the corresponding second variable; The step of each processing circuit weighting the first variable with a first coefficient and adding the first coefficient to the corresponding second variable; The step of each processing circuit using a plurality of the first variables to calculate a problem item; The processing circuits add the problem item to the second variable. The step of each processing circuit calculating a constraint item based on the constraint condition, wherein the constraint item has the first variable as an independent variable; The step of each processing circuit calculating a first correction term including a product of a second coefficient and the constraint term; The step of each processing circuit adding the first correction term to the second variable; The step of each processing circuit increasing the absolute value of the first coefficient and the second coefficient according to the update number; The processing circuits update at least a portion of the first vector and at least a portion of the second vector in parallel; and A step of transferring the first vector between the common storage unit and each of the plurality of processing circuits.

15. A non-transitory computer-readable storage medium, wherein: A program for causing a computer to execute the information processing method according to claim 14 is stored.

16. A computer program product, wherein: A program for causing a computer to execute the information processing method according to claim 14 is included.

Images