Optimization method, optimization device, computer program, and method for solving route optimization problem

Abstract

This optimization method for retrieving a minimum energy state that provides the minimum energy of a spin system includes a step in which a computer retrieves the minimum energy of the spin system when the spin state of a spin included in the spin system is changed. In the retrieval step, the computer divides spins of the spin system into a plurality of spin groups on the basis of causal relationships between the spins, and retrieves the minimum energy of each spin group by changing the spin state of the spins included in the spin group and calculating the energy of the spin system.

Description

Optimization method, optimization device, computer program, and method for solving path optimization problems.

[0001] This disclosure relates to optimization methods, optimization apparatuses, computer programs, and methods for solving path optimization problems. This application claims priority under Japanese application No. 2024-175095, filed on 4 October 2024, incorporating all the provisions of the said Japanese application.

[0002] One type of quantum-inspired optimization method is a technique called quantum annealing. Quantum annealing is specifically designed to solve combinatorial optimization problems. Quantum computers using the quantum annealing method are already partially in practical use. Quantum annealing is a metaheuristic method and does not guarantee the output of the optimal solution. However, compared to other quantum computing methods such as quantum gates, it is expected to output the optimal solution with a higher probability at the current level of technology.

[0003] In quantum annealing, a spin system model called the Ising model is used. This Ising model is described as an undirected graph containing multiple nodes, each assigned a spin, and edges between these nodes. Quantum annealing determines the spin arrangement that gives the minimum energy of the Ising model, assuming that the spins of these nodes are superimposed quantum mechanically. Path optimization problems are a type of combinatorial optimization problem that is difficult to solve, but it is known that candidate solutions can be obtained by determining the spin arrangement in the minimum energy state of such an Ising model.

[0004] One method for simulating this quantum annealing using classical computers (SA (Simulated Annealing)) is called the quantum Monte Carlo method. For the Ising model in quantum annealing, a method called the Suzuki-Trotter decomposition can be used. This method introduces a concept called imaginary time, allowing the Ising model, which is represented by a two-dimensional graph, to be analyzed by a three-dimensional classical model with an increased dimension. Each of these layers is called a Trotter layer. By applying the Monte Carlo method to these Trotter layers, the spin arrangement that gives the minimum energy of the original Ising model can be estimated. Such a technique is disclosed in Non-Patent Document 1, which is listed below.

[0005] Jang Hyung-suk, et al., "Experimental Evaluation of the Path-Integral Monte Carlo Method for Restricted Boltzmann Machines", The 104th Symposium on Fundamental Problems of Artificial Intelligence, The Japanese Society for Artificial Intelligence, Next: 104, Date: 2017 / 08 / 08 - 2017 / 08 / 09, [online], [searched on June 1, 2023], Internet <URL: https: / / www.jstage.jst.go.jp / article / jsaifpai / 104 / 0 / 104_02 / _article / -char / ja / >Kenichi Masuda, et al., "Optimization of Distribution Planning by Quantum Computing", Sumitomo Electric Technical Review, No. 202, January 2023, p. 1-4Hasitha Muthumala Widyasooriyam, et al., "A GPU-based Quantum Annealing Simulator for Fully-connected Ising Models Utilizing Spatial and Temporal Parallelism", [online], IEEE Access, Vol. 8, 9057502, 2020, p. 67929-67939. Internet <url: https: / / doi.org / 10.1109 / ACCESS.2020.<2985699>Hitoshi Goto, "Principles and Applications of Simulated Bifurcation Machines", Surfaces and Vacuum, The Japan Society of Surface and Vacuum Science, Vol. 63, No. 3, 2020, p. 129-133K. Tatsumura, A. R. Dixon and H. Goto, "FPGA-Based Simulated Bifurcation Machine", 2019 29th International Conference on Field Programmable Logic and Applications (FPL), IEEE, 2019, pp. 59-66, doi: 10.1109 / FPL.2019.00019.

[0006] An optimization method relating to the first aspect of this disclosure is an optimization method for searching for a minimum energy state that gives the minimum energy of a spin system, comprising the step of a computer searching for the minimum energy of the spin system when the spin state of the spins included in the spin system is changed, wherein the computer searches for the minimum energy by dividing the spins of the spin system into a plurality of spin groups based on the causal relationships between the spins, and for each of the spin groups, calculating the energy of the spin system by changing the spin state of the spins included in the spin group.

[0007] Figure 1 is a schematic diagram of a two-dimensional Ising model. Figure 2 is a conceptual diagram showing the state in the imaginary time direction after performing a Suzuki-Trotter decomposition on a two-dimensional Ising model. Figure 3 is a flowchart showing the structure of the main routine of a program that implements a path optimization problem solving device employing a quantum inspired optimization method according to one embodiment of this disclosure. Figure 4 is a schematic diagram showing a graph illustrating an example of a causal relationship between spins. Figure 5 is a graph showing an equation representing the relationship between spin vectors when the causal relationship between spins is represented as a Hamiltonian matrix. Figure 6 is a schematic diagram showing a connection matrix when the causal relationship between spins is represented graphically. Figure 7 is a schematic diagram showing in graphical form an example of a spin in the first layer that is directly affected by a certain spin. Figure 8 is a flowchart showing the control structure of a program that implements the quantum inspired optimization of the program shown in Figure 3. Figure 9 is a graph showing the control of temperature and transverse magnetic field in the quantum inspired optimization method according to the first embodiment of this disclosure. Figure 10 is a flowchart showing the control structure of a program that implements the steps of updating the minimum energy and the spin state of the Ising model in the program shown in Figure 8. Figure 11 is a graph showing the results of solving the MaxCut problem using the method according to the first embodiment of this disclosure, along with the results of a comparative technology. Figure 12 is a flowchart showing a part of the control structure of a program that implements a modified version of the first embodiment. Figure 13 is a schematic diagram showing an example of spins up to the second layer that are affected by a certain spin, in graph form. Figure 14 is a flowchart showing the control structure of a program that implements the quantum inspired optimization method according to the second embodiment of this disclosure. Figure 15 is a graph showing the results of solving the MaxCut problem using the method according to the second embodiment of this disclosure, along with the results of the first embodiment and the results of a comparative technology. Figure 16 is a graph showing the degree of speedup due to parallelization when solving the MaxCut problem using the methods according to the first and second embodiments. Figure 17 is a flowchart showing the control structure of a program that implements the quantum inspired optimization method according to the third embodiment of this disclosure.Figure 18 is a flowchart showing the control structure of a program that implements the quantum-inspired optimization method according to the fourth embodiment of this disclosure. Figure 19 is a flowchart showing the control structure of a program that implements the quantum-inspired optimization method according to the fifth embodiment of this disclosure. Figure 20 is a graph showing the degree of speedup of the quantum-inspired optimization method according to this disclosure, compared to when only the Trotter layer is parallelized. Figure 21 is a diagram showing the degree of speedup of the quantum-inspired optimization method according to the first, third, and fourth embodiments of this disclosure, in tabular form for each best cut achievement rate relative to BKC (Best Known Cut), with the case where only the Trotter layer is parallelized as the baseline. Figure 22 is a graph showing the information shown in Figure 21, with the best cut achievement rate relative to BKC on the horizontal axis and the degree of speedup due to parallelization on the vertical axis. Figure 23 is a graph showing the degree of speedup of the first embodiment and the embodiment employing the nested structure shown in the third embodiment in the first embodiment of this disclosure, with the case where only the Trotter layer is parallelized as the baseline. Figure 24 is a table showing the degree of speed improvement for each BKC attainment rate, based on the first embodiment of this disclosure and the embodiment employing the nested structure shown in the third embodiment in the first embodiment, with the BKC attainment rate being the baseline when only the Trotter layer is parallelized. Figure 25 is a graph showing the information shown in Figure 24, with the BKC attainment rate on the horizontal axis and the degree of speed improvement due to parallelization on the vertical axis. Figure 26 is a graph showing the degree of convergence of the optimization calculation when the number of Trotter layers is changed while keeping the number of calculations within one Monte Carlo step (MCS (Monte Carlo Step)) constant, in the third and fourth embodiments, with the MCS number on the horizontal axis and the BKC attainment rate on the vertical axis. Figure 27 is a graph showing the portion of the graph shown in Figure 26 up to MCS 5000, with the MCS value obtained by dividing by the number of Trotter layers used in each graph on the horizontal axis. Figure 28 shows the appearance of a computer operating as a quantum-inspired optimization apparatus, programmed to perform quantum-inspired optimization methods according to each embodiment of this disclosure.Figure 29 is a block diagram showing the hardware configuration of the computer system whose appearance is shown in Figure 28.

[0008] [Problems this disclosure aims to solve] However, conventional quantum Monte Carlo calculation methods have problems such as requiring long computation times and having poor convergence.

[0009] Therefore, this disclosure aims to provide an optimization method, an optimization apparatus, a computer program, and a method for solving a path optimization problem that can reduce computation time and improve convergence.

[0010] [Effects of this disclosure] This disclosure provides an optimization method, an optimization apparatus, a computer program, and a method for solving a path optimization problem that can reduce computation time and improve convergence.

[0011] [Description of Embodiments of This Disclosure] Embodiments of this disclosure are listed below. In the following description and drawings, the same parts are given the same reference numerals. Therefore, a detailed description of them will not be repeated. At least some of the embodiments described below may be arbitrarily combined.

[0012] (1) An optimization method relating to the first aspect of this disclosure is an optimization method for searching for a minimum energy state that gives the minimum energy of a spin system, comprising the step of a computer searching for the minimum energy of the spin system when the spin state of the spins included in the spin system is changed, wherein the computer searches for the minimum energy by dividing the spins of the spin system into a plurality of spin groups based on the causal relationship between the spins, and for each of the spin groups, by changing the spin state of the spins included in the spin group and calculating the energy of the spin system.

[0013] This configuration reduces computation time in the optimization method and improves convergence.

[0014] (2) In (1) above, the search step may include: a step in which the computer identifies a spin group including a spin and a spin affected by the change in the spin state of the spin when the spin state of the spin system is changed by the computer; a state update step in which the computer selectively performs, for each of the spin groups, a process of updating the minimum energy state with respect to the spins included in the spin group, according to whether the energy of the spin system when the spin state of the spin is changed is less than the minimum energy value, and a process of not updating the minimum energy state; an environment update step in which the computer updates the environment of the spin system in response to the completion of the state update step; and an iterative execution step in which the computer repeatedly performs the first process including the step of identifying the spin group, the state update step and the environment update step until a predetermined termination condition is met.

[0015] This configuration reduces computation time and improves convergence in quantum annealing.

[0016] (3) In (2) above, the step of identifying the spin group may include a selection step in which the computer selects one or more arbitrary spins of the spin system as target spins; a step in which, for each of the target spins selected in the selection step, the computer determines which spins in the spin system are directly affected by the target spin; and a step in which the computer identifies the group including the target spins and the spins determined in the determination step as the spin group.

[0017] This configuration allows for parallelization and offline processing of iterative execution.

[0018] (4) In (2) or (3) above, the determining step may include a determination step in which the computer determines the spins that are directly affected by the target spin from the matrix elements of the Hamiltonian matrix describing the spin system that have non-zero values.

[0019] With this configuration, the iterative execution process for spins belonging to non-overlapping spin loops can be parallelized.

[0020] (5) In the above (4), the determining step may further include a step in which a computer uses the spin determined in the determining step as the target spin in the spin system and obtains spins that are directly affected by the target spin from matrix elements having non-zero values in the Hamiltonian matrix.

[0021] With this configuration, the calculation time is further shortened and the convergence is also further improved.

[0022] (6) In any one of the above (2) to the above (5), the first process may include a step in which a computer executes a predetermined second process a predetermined number of times. The second process includes a first iterative step in which a computer repeatedly executes the step of identifying the spin loop and the state update step and ends in response to the establishment of a predetermined first end condition, and a second iterative step in which a computer repeatedly executes the first iterative step and the environment update step and ends in response to the establishment of the end condition.

[0023] With this configuration, the convergence can be enhanced.

[0024] (7) In the above (6), the first iterative step may include a step in which a computer repeatedly executes the step of identifying the spin loop and the state update step and ends in response to the number of repetitions reaching a predetermined number.

[0025] With this configuration, the solution to the energy problem of the spin system at the time when the calculation results have converged to a certain extent can be obtained early.

[0026] (8) In the above (6) or the above (7), the second iterative step may include a step in which a computer repeatedly executes the first iterative step and the environment update step and ends in response to the number of repetitions of the first iterative step and the environment update step reaching a predetermined number.

[0027] With this configuration, the convergence of the quantum-inspired optimization method can be enhanced.

[0028] (9) In any one of (2) to (8) above, the spin system includes a plurality of Trotter layers, and the state update step is such that the computer, in the plurality of Trotter layers, determines whether the energy of the spin system when the spin state of each spin included in the spin group is changed is smaller than the value of the minimum energy calculated by the plurality of Trotter layers, and selectively executes a process of updating the minimum energy state and a process of not updating the minimum energy state. The Trotter layer-by-layer state update step includes this determination.

[0029] With this configuration, the calculation speed of the optimization method can be increased and the convergence can be enhanced.

[0030] (10) In any one of (2) to (8) above, the spin system includes one Trotter layer, and the state update step may include a step in which the computer, in the one Trotter layer, selectively executes a process of updating the minimum energy state and a process of not updating the minimum energy state according to whether the energy of the spin system when the spin state of each spin included in the spin group is changed is smaller than the value of the minimum energy calculated by the one Trotter layer.

[0031] With this configuration, optimization can be performed at high speed using fewer computing resources.

[0032] (11) In (6) above, the first iteration step may include a step in which the computer identifies the spin group in each Trotter layer of the spin system, and a step in which the computer executes the state update step for each spin included in the spin group in all Trotter layers in the step of identifying the spin group.

[0033] With this configuration, the calculation speed of the optimization method can be increased and the convergence can be enhanced.

[0034] (12) In any one of (2) to (5) above, the spin system comprises one or more Trotter layers, and the first process may include a second iteration step in which the computer performs a predetermined second process a predetermined number of times, the second process may include a third iteration step in which the computer performs a predetermined third process a number of times corresponding to the number of the one or more Trotter layers, and a step in which the computer performs the environment update step after the completion of the third iteration step, the third process may include a fourth iteration step in which the computer repeatedly performs the step of identifying the spin group and the state update step and terminates in response to a predetermined first termination condition being met, and a fifth iteration step in which the computer repeatedly performs the fourth iteration step and the environment update step and terminates in response to a termination condition being met.

[0035] This configuration allows for faster computation of the optimization method and improved convergence.

[0036] (13) In (12) above, the state update step may include a step in which the computer selectively performs a process of updating the minimum energy state each time the energy of the spin system when the spin state of each spin included in the spin group is changed in the one or more Trotter layers is less than the minimum energy value calculated by the one or more Trotter layers, and a process of not updating the minimum energy state.

[0037] With this configuration, when the calculated minimum energy value becomes smaller than the previous value, the minimum energy state is updated at that point. As a result, the calculation speed can be increased and convergence can be improved.

[0038] (14) In the above (12) or (13), the number of Trotter layers provided in the spin system may be 1.

[0039] This configuration allows for high-speed optimization using minimal computing resources.

[0040] (15) In the above (12) or (13), the number of Trotter layers provided in the spin system may be multiple.

[0041] This configuration allows for faster computation and improved convergence through parallelization when abundant computing resources are available.

[0042] (16) A computer program relating to the second aspect of this disclosure causes a computer to execute the optimization method described in any one of (1) to (15) above.

[0043] This configuration allows computers to solve combinatorial optimization problems at high speed.

[0044] (17) An optimization apparatus relating to the third aspect of this disclosure includes a computer programmed to perform the optimization method described in any one of (1) to (15) above.

[0045] This configuration allows the computer to be used as a quantum-inspired optimization device capable of solving combinatorial optimization problems at high speed.

[0046] (18) A method for solving a path optimization problem relating to the fourth aspect of this disclosure includes the steps of: a computer converting a path optimization problem formulated in a predetermined format into an energy problem of a spin system; a computer determining a solution for the spin system by applying the optimization method described in any one of (1) to (15) above to the spin arrangement of the spin system; and a computer converting the solution into a solution for the path optimization problem.

[0047] This configuration reduces the computation time required for optimization and improves convergence.

[0048] [Details of Embodiments of this Disclosure] Specific examples of optimization methods, optimization apparatuses, computer programs, and methods for solving path optimization problems according to embodiments of this disclosure are described below with reference to the drawings. This disclosure is not limited to these examples, but is provided for by the claims, and all modifications within the meaning and scope of the claims are intended to be included.

[0049] Section 1 Underlying Technology 1.1 Referring to Figure 1 of the Ising Model, the two-dimensional Ising model 50 includes, as mentioned above, multiple nodes 60, 62, 64, etc., and multiple edges connecting these nodes. Each of the nodes 60, 62, 64, etc., is assigned a spin 70, spin 72, etc. In classical physics, spin is one of the values ​​{1, -1}, representing upward and downward, respectively. On the other hand, in quantum mechanics, it is treated as a superposition of both.

[0050] Figure 2 shows the Ising model 50 after Suzuki-Trotter decomposition. Referring to Figure 2, the Ising model 90 after Suzuki-Trotter decomposition consists of multiple Trotter layers 100, 102, ..., 104. These Trotter layers 100, 102, ..., 104 are spin arrangements decomposed along an axis called the Trotter axis (also called the imaginary time axis). In conventional quantum Monte Carlo methods, the sum of the energies of these Trotter layers 100, 102, ..., 104 was calculated using a classical computer to find the minimum energy state.

[0051] First, the Hamiltonian H representing the energy of the Ising model 50 system is expressed by the following equation:

[0052] In this formula, J ｉｊ This indicates the magnitude of the interaction between adjacent spins in the spin array set as the problem. (σ with hat) ｉ ｘ and σ ｉ ｚ These are Pauli spin matrices, and are represented by the following matrices:

[0053] Γ represents the magnitude of the transverse magnetic field applied to the Ising model.

[0054] The partition function Z can be obtained from this Hamiltonian H by the following equation (3). However, in equation (3), generally β = 1 / kＢ T(k Ｂ is the Boltzmann constant, and T is the temperature). |ψ〉 is a state vector on the Hilbert space representing the state of the system. In SA, β varies in proportion to the reciprocal of the temperature. Therefore, any control is allowed in the calculation. In the first embodiment, this β is made constant.

[0055] The goal is to find the energy E that maximizes this partition function Z, that is, to find the spin arrangement that minimizes the energy.

[0056] When the Suzuki-Trotter decomposition is applied to the above formula (1), it can be represented by a three-dimensional classical representation in which the imaginary time axis is introduced.

[0057] In this formula, P represents the number of Trotter layers, and s ｋ ｉ represents the i-th spin (1 or -1) of the k-th Trotter layer. J ⊥ represents the magnitude of the interaction between spins in different (adjacent) Trotter layers, and T represents the temperature in annealing.

[0058] As a result, as described above, solving the combinatorial optimization problem reduces to finding the state that gives the minimum value of the energy E in the Ising model, that is, the spin arrangement.

[0059] 1.2 Calculation Method In ordinary quantum-inspired optimization that faithfully calculates this formula, as represented by the following formula, it is necessary to calculate the sum of the energies of each Trotter layer.

[0060] When calculating the above formula, in the prior art, it was necessary to use the following method to obtain the minimum energy value. That is, the spin arrangements are updated independently for all of the Trotter layers 100, Trotter layer 102,..., Trotter layer 104 shown in FIG. 2. Further, using the result of adding them together, the spin state and the minimum value of the energy of the Ising model 50 are updated.

[0061] However, this calculation requires performing the calculation shown in equation (7) for all Trotter layers. Therefore, in conventional techniques, there was a problem in that the time taken for each iteration of the iterative calculation became long and convergence was poor.

[0062] To solve these problems, the applicant of this patent has already proposed the following: perform the calculation of argmin first, and omit the process of calculating the sum of all energies in the Trotter layer.

[0063] In this equation, |α Ｋ > is a new spin sequence candidate. This new spin sequence candidate is the K-th spin sequence |α obtained immediately before. Ｋ In >', a predetermined number of spins (one in this embodiment) are randomly selected and inverted to obtain the result. On the other hand, |α ｍｉｎ > is the spin arrangement that gives the minimum energy at the current time out of all P spin arrangements in the Trotter layer, and is a hat-shaped σ ｚ ｍｉｎ，ｉ E1 is the Pauli matrix acting on the spin array that gives the minimum energy. Therefore, E1 represents the energy given by the interaction between spins in the K-th Trotter layer. E2 represents the energy given by the interaction between the spin array of the K-th Trotter layer and the spin array that gives the minimum energy at the present moment.

[0064] According to this formula, the energy is not calculated for the entire Trotter layer before finding the minimum value. In this formula, the spin arrangement |α that minimizes the energy in the Trotter layer is... ｍｉｎ If a minimum value is found, that value is immediately used in the calculation of other spin arrangements. In this way, calculations for each Trotter layer can be performed in parallel, and when the minimum value is found, that value can be immediately reflected in the calculations for each Trotter layer. Therefore, the calculation method shown in this embodiment reduces the amount of computation required to obtain the spin arrangement that gives the minimum energy. As a result, the solution can be obtained at high speed. The convergence to the spin arrangement that gives the minimum energy is also improved.

[0065] This disclosure builds upon and further improves upon the methods described above.

[0066] 2.1 Configuration The quantum-inspired optimization according to this first embodiment can be applied to various fields. One example is the optimization of delivery plans for packages in the logistics field. Delivery plan optimization is a combinatorial optimization problem, but the computational complexity increases exponentially as the number of vehicles used and the number of delivery destinations increases, making it difficult to find an exact solution. The quantum-inspired optimization described above has the potential to significantly reduce the time required to solve the delivery plan optimization problem.

[0067] As described in Non-Patent Document 2, the above-mentioned problem concerning the optimization of delivery plans can be formulated using a method called QUBO (Quadratic Unconstrained Binary Optimization). Furthermore, it is known that QUBO can be transformed into a graph max-cut problem. The max-cut problem is the problem of finding a way to divide a graph into two subgraphs that maximizes the number of edges that are cut (maximization is one example of optimization). The max-cut problem is known to be equivalent to the Ising model. Therefore, the quantum-inspired optimization method described above can be applied to the max-cut problem.

[0068] Figure 3 is a flowchart showing the control structure of the main routine of a delivery plan optimization program to which the quantum-inspired optimization method according to this first embodiment is applied. Referring to Figure 3, the program includes a step 150 of reading a QUBO formulation file that records a delivery plan optimization problem formulated in QUBO format, and a step 152 of converting the contents of the formulation file read in step 150 into a max-cut problem and an Ising model equivalent to this max-cut problem. The program further includes a step 154 ​​of determining the spin arrangement of the minimum energy state of the Ising model by applying the quantum-inspired optimization method according to this disclosure to the Ising model obtained in step 152, a spin update table 158 that tables the causal relationships between spins and is referenced in the processing in step 154, and a step 156 of outputting the spin arrangement calculated in step 154 ​​as the answer to the delivery plan optimization problem formulated as the initial QUBO problem and ending the processing.

[0069] The causal relationships between spins will be explained with reference to Figures 4 and 5. Referring to Figure 4, consider a spin group 190 consisting of multiple spins. When the state of one spin in spin group 190 changes, the state of one or more spins within spin group 190 changes as a result. This can be considered as a causal relationship existing between each spin included in spin group 190 and the other spins in the same spin group 190. This causal relationship is shown in graph form in Graph 180 shown in Figure 4.

[0070] In other words, graph 180 includes spin group 190, spin group 192 consisting of the same elements as spin group 190, and a set of directed edges (edge ​​set 194) from each spin of spin group 190 to one or more spins of spin group 192. Each element of spin group 190 represents a causal relationship from the state of a spin that is an element of spin group 190 to the state of other spins that are elements of spin group 192.

[0071] The relationship shown in Graph 180 represents an example of a Markov process, as it indicates that a new state of spin group 190, spin group 192, is determined based on spin group 190.

[0072] This relationship is clearly reflected in the shape of the Hamiltonian. Referring to Figure 5, consider a vector 222 whose elements are values ​​representing the states of each element of the spin group, and a vector 220 related to vector 222 by the Hamiltonian 224. When only the j-th element vj of vector 222 is changed, the Hamiltonian 224 determines which element of vector 220 is affected.

[0073] Causal relationships can also be represented by the form of the connection matrix shown in graph 180 in Figure 4. Figure 6 shows an example of such a connection matrix, connection matrix 226. Referring to Figure 6, let N be the number of elements included in spin group 190. The connection matrix 226 is an N x N square matrix. The value of the cell in row i, column j of the connection matrix 226 indicates whether or not a causal relationship exists from the i-th spin to the j-th spin. From this matrix, the spin update table 158 shown in Figure 3 is created and stored in the memory.

[0074] Figure 7 schematically shows, in graph form, a spin group 230 whose elements are a spin 240 and spins that are affected by the state change of that spin (i.e., have a causal relationship with spin 240). In this embodiment, a spin group that is directly affected by the state change of a spin, such as spin group 230, is called the "first layer". The connection matrix 226 shown in Figure 6 stores information for all spins such that, for example, given a spin 240, the spin group 230 in the first layer that is affected by the state change of spin 240 can be identified.

[0075] Referring to Figure 8, step 154 ​​shown in Figure 3 involves performing a Suzuki-Trotter decomposition on the Ising model corresponding to the Maxcut problem, thereby obtaining an initial spin arrangement (|α) of a predetermined number (let's call it P) of spin arrangements. ＫStep 300 involves randomly generating and preparing each of the following: |α (where K=1 to P), and setting the initial state of the spin arrangement that gives the minimum energy of the entire system as the |α Ｋ The process includes step 302 of searching for a spin arrangement that gives the minimum energy, and step 304 of setting initial values ​​for temperature (T) and transverse magnetic field strength (Γ) to be used as conditions for energy calculations in the following iterations.

[0076] Step 154 ​​further includes step 306, which repeats step 308, a process of calculating the minimum energy for all trotter layers while keeping the temperature T used as a condition for energy calculation constant and changing the transverse magnetic field Γ, for a predetermined number of times (i.e., until the termination condition is met), and step 310, which, at the end of step 306, inversely transforms the spin arrangement that gives the minimum energy of the Ising model obtained during the process of step 306 into the solution to the delivery plan optimization problem formulated as QUBO and outputs it, thereby terminating the process.

[0077] In this first embodiment, the temperature T is kept constant, as shown by the straight line 384 in Figure 9. The transverse magnetic field Γ is updated, for example, with each execution of step 308, from a predetermined initial value to a very small value after a predetermined number of iterations as described above, as shown by the curve 380. In this embodiment, the value of the transverse magnetic field Γ is controlled to be inversely proportional to the number of Monte Carlo steps (MCS). The J obtained by equation (4) ⊥ The value of is shown by curve 382. This shows that J ⊥ The graph is obtained as curve 382 shown in Figure 9. As a result, as shown by curve 382, ​​J is responsible for the effect of aligning the spin arrangement of the Trotter layer. ⊥ The contribution of this becomes larger.

[0078] Step 308, shown in Figure 8, includes step 330, which repeats step 332, which updates the minimum energy and spin arrangement, for the number of trotter layers while changing the spin arrangement of the target trotter layer, and step 334, which, after the completion of step 330, updates the values ​​of temperature (T) and transverse magnetic field (Γ) according to Figure 9, and then terminates the execution of step 308.

[0079] Step 332 includes step 350, which randomly selects a number of spins from the spin arrangement of the target Trotter layer. The number of spins selected in step 350 may be one. Step 332 further includes step 354, which performs energy calculations of the Ising model and updates the state of the Ising model while changing the state of each spin in a spin group that includes each spin selected in step 350 and the spins of the first layer that are related to that spin in the spin group 230 shown in Figure 7.

[0080] Step 354 includes step 370, which determines the spin in the first layer for each of the spins selected in step 350 by referring to the spin update table 158 shown in Figure 3, and determines a spin group together with the spins selected in step 350, and step 372, which performs an energy calculation state update process 374 for each spin included in this spin group.

[0081] Figure 10 shows a more detailed control structure of the state update process 374 in Figure 8 in flowchart format. Referring to Figure 10, the state update process 374 includes a step 400 to generate a new state candidate by changing the state of the spin to be processed, a step 402 to calculate the Ising model energies E1 and E2 for this state candidate, and a step 404 to calculate ΔE1 and ΔE2, which are the differences between the new energies E1 and E2 calculated in step 402 and the energies E1 and E2 calculated from the current spin arrangement.

[0082] The state update process 374 further includes a step 406 in which the control flow is branched according to whether ΔE1 is 0 or less, and in response to the determination in step 406 being positive, the current spin arrangement |α Ｋ >' is a candidate state for the spin arrangement |α ＫThe state update process 374 further includes step 408, which updates the state to >, and step 410, which follows step 408, which branches the control flow according to whether the energy E1 calculated in step 402 is less than or equal to the minimum energy Emin that has been updated by previous calculations. In this specification, the terms "less than or equal to" and "less than or equal to" may be interpreted interchangeably. The state update process 374 further updates the minimum energy Emin with the new energy E1 in response to the determination in step 410 being positive, and sets the spin array that gives the minimum energy Emin to the spin array state candidate |α Ｋ The process includes step 412, which updates the state and terminates the state update process 374. If the determination in step 410 is negative, the execution of the state update process 374 is terminated.

[0083] The state update process 374 further includes step 414, in response to a negative determination in step 406, a step 414 that branches the control flow according to whether ΔE1 + ΔE2 (= difference of E1 + difference of E2) is 0 or less; step 416, in response to a positive determination in step 414, a step 416 that samples a random number ξ from a uniform distribution of [0,1]; step 418, in response to a positive determination in step 416, a step 418 that branches the control flow according to whether the value ξ sampled in step 416 is less than or equal to exp(-λ(ΔE1 + ΔE2) / T); and in response to a positive determination in step 418, the spin array |α to be processed. Ｋ >' is a candidate state | α Ｋ The procedure includes step 420, which updates the state and terminates the execution of the state update process 374. If the determination in step 414 is negative, and if the determination in step 418 is negative, the state change is not performed and the execution of the state update process 374 is terminated.

[0084] In this embodiment, in step 406, the state change is determined based solely on ΔE1. ​​This is because all problems to be solved are set to Jij, and therefore the focus is on minimizing the energy E1 that includes Jij. Here, Ji is the i-th row vector of the coupling coefficient matrix J between spins, and Jij is its j-th element. However, this disclosure is not limited to such embodiments. For example, in step 406, the control flow may be branched depending on whether the value of the linear sum of ΔE1 and ΔE2 (a × ΔE1 + b × ΔE2) is less than or equal to 0. The above-described step 406 is an example where a = 1 and b = 0 in this linear sum.

[0085] Thus, the coefficients a and b in the above linear combination are determined by appropriately reflecting the contribution of ΔE² in the above calculation. In other words, these values ​​can also be determined such that a + b = 1. In this case, it is not the absolute values ​​of coefficients a and b that are significant, but the ratio between them. Thus, the value of (a, b) can be a = b = 1, a + b > 1, or a + b < 1. Furthermore, a × ΔE １ ｎ +b × ΔE² ｍ This form is also acceptable, and can be arbitrarily selected depending on the Ising model to be solved.

[0086] In this way, even if the energy E1 increases due to the candidate state, if it approaches the spin arrangement that minimizes the energy, the spin arrangement of the Trotter layer is replaced with the candidate state with a certain probability. As a result, the probability that the energy will remain in the local minimum can be reduced while maintaining convergence. The λ that appears in the above equation can be adjusted as appropriate to parameters such as the number of Trotter layers. In addition, the routine shown in Figure 10 can be modified in various ways. For example, step 408 may be moved to a location where it is executed when the determination in step 410 is affirmative. In step 410, only ΔE1 is considered, but ΔE1 + ΔE2 may be considered. In step 418, ΔE1 + ΔE2 is considered, but ΔE1 only may be considered. As mentioned above, in the comparison step, "less than or equal to" may be replaced with "less than". The parameters used in this routine can also be flexibly changed within a range that does not cause theoretical inconvenience.

[0087] At this point, if the minimum energy and the spin array that provides the minimum energy are updated in step 412, the spin array that provides the minimum energy may be used in step 412 of the processing of other Trotter layers in step 332. There is no need to wait for the end of step 332. Step 330 executes step 332 repeatedly. In this embodiment, as described above, the control structure is such that when a spin array that provides the minimum energy is found, the spin array that provides the minimum energy is updated without waiting for the end of step 332. Therefore, if the iterative processing of step 332 is performed in parallel and a spin array that provides the minimum energy is obtained in any of the processes, that spin array can be used in all other processes to perform subsequent processing. In other words, this iterative processing can be executed in parallel on a GPU (Graphics Processing Unit) or the like, and the time until a solution is obtained can be shortened.

[0088] In steps 402, 404 and 406 to 420, for each of a predetermined number of spin arrays, the state of the trotter layer and the minimum energy state are updated based on the energy of the state candidate, the energy of the spin array to be processed, and the minimum energy.

[0089] As described later, the quantum-inspired optimization method described above requires less time and has a higher probability of obtaining the optimal solution in the benchmark problem of the max-cut problem compared to conventional methods. This is because, in the state update process for spin arrangements, if there is a spin arrangement that gives the minimum energy, that value is immediately used in subsequent calculations. As a result, the computational load is reduced because the number of iterative loops is effectively reduced by one compared to the conventional method of selecting the spin arrangement that gives the minimum energy after calculations for all spin arrangements are completed.

[0090] Furthermore, in this embodiment, when one spin is selected, the spin arrangement can be explored not only for that spin but also for the spins in the first layer that are affected by the state change of that spin. As a result, as described in the benchmark section below, simulations confirmed that the processing time is significantly reduced compared to the case where the spins in the first layer are not used, that is, when the state update process 374 is performed only on spins randomly selected from the same number of spins instead of step 370 in Figure 8.

[0091] 2.2 Referring to Operation Diagrams 8 and 10, the computer for implementing this quantum-inspired optimization method operates as follows. This computer operates as a quantum-inspired optimization device by being programmed with the programs shown in Figures 8 and 10.

[0092] Referring to Figure 8, first, in step 300, by performing a Suzuki-Trotter decomposition on the Ising model corresponding to the Maxcut problem, the arrangement of P Trotter layers |α Ｋ>' is prepared randomly. In step 302, the initial state of the spin array that gives the minimum energy of the entire system is |α Ｋ >' is searched for. Any method can be used for this search. |α Ｋ The value of 〉' may be set to a predetermined initial value. Subsequently, in step 304, the temperature (T) and the strength of the transverse magnetic field (Γ), which are used as conditions for energy calculation in the following iterations, are set to predetermined initial values.

[0093] Subsequently, in step 306, step 308 is repeated a predetermined number of times. Since this process is performed on different Trotter layers, it can be executed in parallel. In this iterative process, the following processes are performed. First, step 332 is performed for each Trotter layer. That is, step 332 is performed for each of the different Trotter layers.

[0094] In step 350 of step 332, multiple spins within the Trotter layer to be processed are randomly selected. However, as mentioned above, only one spin may be selected in step 350.

[0095] In the subsequent step 352, step 354 is performed for each of the multiple spins selected in step 350. That is, step 354 is performed as many times as there are spins selected in step 350. Since this process is performed for each selected spin, it can be executed in parallel.

[0096] Inside step 354, first, the spin group affected by the change in the state of the selected spin, i.e., the spin group of the first layer, is determined. For this determination, the elements of Hamiltonian 224, shown in the equation illustrating the relationship between spins in Figure 5, are used. That is, the elements of vector 222 corresponding to the selected spin are related to each element of vector 220 by Hamiltonian 224. Therefore, by examining the shape of each element in the column of Hamiltonian 224 that corresponds to the element of vector 222 corresponding to the selected spin (e.g., vj), it can be determined which elements of vector 220 are affected by which elements of vector 22. In this way, for example, the spin of the first layer affected by the spin corresponding to element vj is determined as the spin in vector 220 that is affected by the change in element vj.

[0097] The spins of the first layer determined in step 370 and the spins being processed in step 354 (one of the spins selected in step 350) can be considered as a single group. In step 372, the state update process 374 is performed for each spin belonging to this group.

[0098] Figure 10 shows the details of the state update process 374. Referring to Figure 10, the following processes are performed in the state update process 374. In step 400, the spin of the material to be processed is reversed. This spin determination changes the state of the trotter layer being processed. In step 402, the energy E1 and energy E2 after this state change are calculated. In step 404, the difference ΔE1 and difference ΔE2 are calculated, which are the difference from the energy E1 and energy E2 obtained in the previous execution of the state update process 374.

[0099] In the following step 406, it is determined whether the value of the difference ΔE1 is 0 or less, and the control flow branches according to the result. This determination may also be based on whether the value of the difference ΔE1 is less than 0. If the determination in step 406 is positive, it means that a new minimum energy value has been found. Therefore, the control proceeds to step 408, in which the current spin arrangement |α Ｋ >' is a new state candidate obtained by the processing in step 400 |α Ｋ > will be updated.

[0100] Subsequently, in step 410, the control flow is branched depending on whether the energy E1 is less than or equal to the minimum energy Emin. If the determination in step 410 is positive, in step 412 the minimum energy Emin is updated with the new minimum energy E, and the spin array |α gives the minimum energy Emin. ｍｉｎ > Candidate spin arrangement state | α Ｋ The execution of this step ends with the update. The execution of this step also ends if the judgment in step 410 is negative or positive.

[0101] On the other hand, if the determination in step 408 is negative, the control proceeds to step 414. In step 414, the control flow is branched depending on whether ΔE1 + ΔE2 is less than or equal to 0. If the determination in step 414 is negative, the execution of this step ends. If the determination in step 414 is positive, a random number ξ is sampled from a uniform distribution in [0,1] (step 416), and the control flow is branched depending on whether this value ξ is less than or equal to exp(-λ(ΔE1 + ΔE2) / T) (step 418). If the determination in step 418 is positive, the spin array to be processed |α Ｋ >' is a candidate state | α Ｋ The state is updated (step 420) and the execution of the state update process 374 is terminated. If the determination in step 414 or step 418 is negative, this state change is not performed and the execution of the state update process 374 is terminated.

[0102] Referring again to Figure 8, the state update process 374 is repeated for all spins in the spin group determined in step 370 until it is completed. Once the state update process 374 has been executed for all spins in the spin group, the processing in step 354 for the spin selected as the target of processing in step 352 is completed, the next spin is selected, and the state update process 374 is started again.

[0103] Thus, once step 354 is completed for all the spins selected in step 350, the process in step 332 is performed on the next trotter layer selected in step 330. Furthermore, once the process in step 332 is completed for all trotter layers, the magnetic field and temperature are updated in step 334 according to the schedule, and one execution of step 308 is completed. After repeating the above step 308 a predetermined number of times, step 310 is executed. In step 310, the state that gives the minimum energy among the energies obtained in the above iterative process is selected as the optimal solution.

[0104] The optimal solution obtained in this way is then transformed back into QUBO's answer in step 156 of Figure 3. As a result, by transforming it back into the answer to the initially given delivery plan optimization problem, one answer to the delivery plan is obtained.

[0105] 2.3 Benchmark The results of a benchmark experiment comparing the quantum inspired optimization method according to this embodiment with a quantum inspired optimization method that does not use the spin of the first layer are shown below. The quantum inspired optimization method that does not use the spin of the first layer corresponds to the one already proposed by the applicant of this application, as described in "Underlying Technology".

[0106] This experiment used a dataset called Gset, which is commonly used in this field as a benchmark for the MaxCut problem.

[0107] Figure 11 is a graph showing the results obtained by applying the method based on the prerequisite technology to dataset G34, one of the benchmark datasets of Gset, and the results obtained by applying the method according to the first embodiment described above. In Figure 11, the horizontal axis represents the Monte Carlo steps (MCS), and the vertical axis represents the number of best cuts obtained up to each MCS. Here, in the first embodiment described above, an extra step of several times the spin of the first layer is performed within one MCS. Therefore, the results shown are after correcting for this effect.

[0108] In the graph shown in Figure 11, curve 430 shows the results obtained using the underlying technology, and curve 432 shows the results obtained using the method of the first embodiment. Referring to Figure 11, it can be seen that curve 432 surpasses the results of curve 430 in all MCSs. In other words, this simulation suggests that by using the method of the first embodiment, better results can be obtained in a shorter time compared to the underlying technology.

[0109] 2.4 Modified Examples In the above embodiment, spin group determination is performed internally in step 332. However, this disclosure is not limited to such embodiments. Spin group determination may be performed in advance as an offline process and stored in a memory device as a spin group table. When the optimization method is executed, spins can be read from the spin group table. Figure 12 shows a flowchart of the control structure of a program for implementing such a modified example.

[0110] Referring to Figure 12, the program for creating the spin group table 445 described above offline includes a step 440 of randomly selecting a plurality of spins, and a step 441 of executing step 442 for each of the spins selected in step 440 to store information indicating the spins that constitute the spin group in the spin group table 445.

[0111] Step 442 includes step 443, which determines the spin of the first layer for the spin selected in step 441, and step 444, which writes the spin determined in step 443 and the spin selected in step 441 as a spin group to the spin group table 445. When the repetition of step 441 is completed, the spin group table 445 is stored in a predetermined storage device (not shown).

[0112] The flowchart on the left side of Figure 12 shows a program routine that can be replaced with step 308 of the first embodiment shown in Figure 8.

[0113] Referring to the left side of Figure 12, step 450, which is interchangeable with step 308, includes step 451, which repeats step 452 for a given Trotter layer, for a number of Trotter layers, changing the Trotter layer in question, so as many times as there are Trotter layers, and step 334, which updates the magnetic field and temperature according to a schedule in response to the completion of step 451.

[0114] Step 452 includes step 453 of reading spin group data from spin group table 445, and step 454 of performing a state update process 374 for each spin group read in step 453 to search for a spin configuration that gives the minimum energy of the Ising model.

[0115] As shown in this modified example, even if step 444 is processed offline, the same processing as in Embodiment 1 can be performed. This is also true for other embodiments described later. Offline processing is possible even when the number of spin layers used is increased.

[0116] Third Second Embodiment 3.1 Configuration The first embodiment described above makes it possible to perform optimization in a shorter time than the presupposed technology. However, this disclosure is not limited to the embodiments of the first embodiment. For example, in the first embodiment, a group is formed for a selected spin, including the spins in the first layer that are affected by that spin, and the minimum energy state of the Ising model is searched for these groups. However, the search for the minimum energy state may be performed not only for the first layer, but also for each spin belonging to the first layer, including the spins that have a relationship with the first layer. In other words, from the perspective of the initially selected spin, the minimum energy state can be searched for each spin in the spin group that includes not only the first layer, but also the second layer.

[0117] In this second embodiment, the search for the minimum energy state is performed including the second layer described above.

[0118] Figure 13 graphically shows the relationship between the selected spin 240, the spin group 230 consisting of the first layer spins related to spin 240, and the spin group 460 consisting of the second layer spins related to spin 240. Figure 13 only shows the second layer spin group of spin group 230, consisting of the first layer spins 470, 472, ..., 474 related to spin 242. However, similarly, there are second layer spins for each of the first layer spins 244, ..., 246. Therefore, the number of spins included in spin group 460 is larger than the number of spins included in spin group 230.

[0119] Figure 14 shows, in flowchart form, the control structure of a computer program for causing a computer to implement the quantum-inspired optimization method according to this second embodiment.

[0120] Referring to Figure 14, this program differs from the first embodiment shown in Figures 8 and 10 in that, instead of step 306 in Figure 8, in which step 308 is executed a predetermined number of times, it includes step 480, in which step 482 is executed a predetermined number of times. In all other respects, Figure 14 is the same as Figure 8.

[0121] Step 482 shown in Figure 14 includes step 500, which repeats step 502 as many times as there are Trotter layers, and step 334, which is the same as that shown in Figure 8.

[0122] Step 502 is similar to step 332 shown in Figure 8, but differs in some respects. More specifically, step 502 differs from step 332 in that it includes step 350, which is the same as step 332, and step 510, which replaces step 352 in Figure 8, by performing step 512 for each spin selected in step 350.

[0123] Step 512 includes step 370, which determines the spin of the first layer for the spin to be processed, and step 530, which performs step 532 for each spin determined in step 370.

[0124] Step 532 determines the spins of the first layer (which become the spins of the second layer for the spins selected in step 350) for the spins to be processed from among the spins of the first layer. Upon completion of step 530, a spin group is determined consisting of the spins selected for processing in step 510, the spins of the first layer that were processed in step 370, and the spins of the second layer determined in step 530.

[0125] Step 512 further includes step 372, in which the same state update process 374 as the state update process 374 in Figure 8 is performed for each spin in the spin group determined in step 530.

[0126] In short, this second embodiment has a configuration that differs from the first embodiment, which only performs the state update process 374 on the spins up to the first layer, in that it performs the state update process 374 not only on the spins up to the first layer but also on the spins up to the second layer.

[0127] 3.2 Operation The quantum inspired optimization method according to this second embodiment is substantially the same as that of the first embodiment. The differences are as follows. In the first embodiment shown in Figure 8, a spin group consisting of the spin selected in step 352 and the spin of the first layer is determined by determining the spin of the first layer in step 370 at the beginning of step 354. Then, in step 372, a state update process 374 is performed for each spin in this spin group. In contrast, in the second embodiment shown in Figure 14, step 354 corresponds to step 512. In step 512, as in the first embodiment, the spin of the first layer is selected in step 370 for the spin specified in step 510. In this embodiment, step 532 is further performed for each of the spins of the first layer. In step 532, the spins of the second layer that are affected by the spins of the first layer are determined. Then, in step 372, a state update process 374 is performed for each spin in the spin group consisting of the spin selected in step 510, the spin of the first layer selected in step 530, and each spin selected in step 532. Here, there may be cases where the selected spins overlap, in which case they may be combined into one to avoid duplication.

[0128] In other words, in the first embodiment, for each spin selected in step 350, a spin group is formed including the spins of the first layer relative to that spin. In contrast, in the second embodiment, for each spin selected in step 350, a spin group is formed including the spins of the first layer relative to that spin, and further including the spins that are in a first-layer relationship with the spins of the first layer (the spins of the second layer relative to the spin selected in step 350). This is the major difference between the two.

[0129] To determine the spin of the second layer, the spin update table 158 shown in Figure 3 can be used twice.

[0130] 3.3 Benchmark Figure 15 shows the simulation results for the method based on the underlying technology, the quantum-inspired optimization method according to the first embodiment, and the quantum-inspired optimization method according to the second embodiment. The results for the method based on the underlying technology, the results for the first embodiment, and the results for the second embodiment are shown by curves 580, 582, and 584, respectively. In Figure 15, the horizontal axis represents MCS, and the vertical axis represents the obtained max-cut values. However, this graph does not correct for the number of extra trials included in 1 MCS.

[0131] Referring to Figure 15, this second embodiment approaches BKC at a much smaller MCS time compared to the method of the preceding technology. Compared to the first embodiment, the result shows that BKC is approached much faster. Therefore, when it is sufficient to obtain a solution with a certain degree of accuracy, even if it is not a solution very close to the optimal solution, the quantum-inspired optimization method according to this second embodiment is very effective if the spin selection in step 350 of Figure 14 is parallelized.

[0132] Figure 16 shows the results of the first embodiment and the second embodiment as graphs with the rate of reaching BKC on the horizontal axis and the degree of speed improvement due to parallelization on the vertical axis. The degree of speed improvement referred to here is based on the results when using the comparison technology.

[0133] Referring to Figure 16, in the first embodiment, compared to the baseline, the degree of speed increase increases gradually as the rate of reaching BKC increases. However, even when the rate of reaching BKC approaches 1.0, the rate of speed increase is only about 7 times, and in this example, the average number of spins up to the first layer is about 5. Therefore, the degree of speed increase can be evaluated as being roughly the same as the number of spins selected to be processed at one time. In contrast, in the second embodiment, it can be seen that the degree of speed increase is much higher than in the first embodiment. In particular, the degree of speed increase increases as the rate of reaching BKC approaches 1.0, and at around 1.0, the degree of speed increase reaches 25 times. In this example, the average number of spins up to the second layer processed at one time was 13. Therefore, the degree of speed increase converges faster compared to the first embodiment.

[0134] Compared to the comparable technology, this second embodiment can further improve convergence and enable optimization in a shorter time than the first embodiment.

[0135] 4.3 Third Embodiment 4.1 Configuration Figure 17 shows the control structure of a program for implementing the quantum inspired optimization method according to the third embodiment of this disclosure using a computer. Referring to Figure 17, this program is the same as the program shown in Figure 14, but with the contents of step 482 further modified into a nested structure.

[0136] More specifically, this program includes the same steps 300, 302, and 304 as in Figure 14. This program further includes step 610, which replaces step 480 in Figure 14 by performing step 612 a predetermined number of times, and step 310, which, after the completion of step 610, selects the spin state that gives the minimum energy of the Ising model obtained in step 610 and terminates the process.

[0137] Step 612 includes step 620, which performs the first iteration process 622 a predetermined number of times, and step 334, which, after the completion of step 620, updates the values ​​of temperature (T) and transverse magnetic field (Γ) according to the schedule shown in Figure 9, and then terminates one execution of step 612. This step 334 is the same as step 334 shown in Figures 4 and 14.

[0138] The first iteration 622 includes, as shown in Figures 4 and 14, a step 350 of randomly selecting a plurality of spins, and a step 640 of performing a step 642 for each selected spin to calculate the minimum energy in each trotter layer.

[0139] Step 642 includes, as shown in Figures 4 and 14, a step 370 in which a first layer spin is determined for the spin to be processed; a step 530 in which a spin group to be processed is determined by performing a step 532 in which a first layer spin is determined for each of the first layer spins (the second layer spin from the perspective of the spin selected in step 350); and a step 680 in which a second iterative process 682 is performed in which the energy calculation and state update of the Trotter layer are performed for each of the spins included in the spin group determined in step 350 when the state of that spin is changed.

[0140] The second iteration process 682 includes a step 690 in which the same state update process 374 as shown in Figure 14 is performed for each Trotter layer.

[0141] In other words, in this embodiment, the first process is repeatedly executed in step 620. Within the first iterative process 622, the iterative process in step 640, which selects one spin from among the causally related spins and executes step 642, and the process in the second iterative process 682, which searches for the minimum energy to be performed for each Trotter layer, are nested.

[0142] This third embodiment can be considered a process that focuses on a spin group containing a certain spin in the Ising model and partially optimizes each Trotter layer only for that spin group.

[0143] The effects of this third embodiment will be described later, along with the effects of the fourth embodiment, which will be discussed later.

[0144] 4.2 Operation The operation of this third embodiment is more complex than that of the second embodiment in that it incorporates a nested structure. However, apart from the nested structure, the main aspects of the operation are the same as those of the second embodiment.

[0145] Step 612 in the third embodiment shown in Figure 17 corresponds to step 482 shown in Figure 14. Step 612 includes steps 620 and 334. Step 612 includes the iteration process by step 620, the iteration process by step 642, and the iteration process by step 680. In particular, the nested structure formed by the iteration process by step 620 and the iteration process by step 680 is a feature of this embodiment.

[0146] In this embodiment, after steps 300, 302, and 304 are performed, step 610 performs a predetermined number of repetitions of step 612.

[0147] In the process of step 612, after the first iteration process 622 is executed a predetermined number of times, step 334 is executed.

[0148] In one execution of the first iteration 622, first in step 350, multiple spins are randomly selected. These spins are treated as common to multiple Trotter layers.

[0149] In step 640, step 642 is performed for each of these selected spins. That is, step 642 is performed as many times as there are spins selected in step 350.

[0150] During the execution of step 642, first, in step 370, the spins of the first layer for the spins to be processed are determined. Furthermore, for each of these spins in the first layer, the spins of the second layer are determined in step 532 of step 530. As a result, the spin group to be processed is determined. To obtain the spins of the second layer from the spins of the first layer, the spin update table 158 shown in Figure 3 can be used twice. Steps 370 and 530 determine the spin group consisting of the spins for which energy calculations are to be performed.

[0151] In step 680, the second iteration 682 is performed for each spin belonging to this spin group. In step 690 of the second iteration 682, the state update process 374 is performed for the same spin in each trotter layer. As a result, for each trotter layer, energy calculations and state updates are performed for the state change of the spin being processed in that trotter layer. Once step 690 is completed by performing the state update process 374 for all trotter layers, the same process is performed for other spins in the same spin group by repeating step 680. In this way, once all of the second iterations 682 have been performed for the spin group being processed, one execution of step 642 is completed.

[0152] Next, if there are any other spins remaining among the spins selected in step 350, step 642, consisting of the above processes, is executed in the same manner. Once step 642 has been completed for all selected spins, one iteration of the first iteration process 622 is completed.

[0153] Next, the first iteration 622 is repeated through the iteration process of step 620. In this process, a spin different from the spin selected in the previous iteration is selected in step 350.

[0154] In this way, once the first iteration process 622 has been executed a predetermined number of times, step 620 is completed. Then, in step 334, the values ​​of temperature and transverse magnetic field are updated according to the schedule, and step 612 is executed again.

[0155] Once step 612 has been executed a predetermined number of times, the process in step 610 is terminated. In step 310, the state that provides the minimum energy calculated in the processes up to that point is selected, and the execution of this program ends.

[0156] This third embodiment enables a higher level of parallelization, and computational simulations confirmed that it achieves even greater convergence than the second embodiment. The details will be described later, along with the results of the fourth embodiment. Furthermore, this nested structure, as in the third embodiment, can also be applied when the spin group includes spins up to the first layer.

[0157] Fifth Fourth Embodiment 5.1 Diagram 18 shows, in flowchart form, the control structure of a computer program for causing a computer to implement the quantum inspired optimization method according to the fourth embodiment of this disclosure.

[0158] The fourth embodiment shown in Figure 18 and the third embodiment shown in Figure 17 share a common structure in that the iterative processes are nested. However, the difference between the two lies in the fact that the nested structure, consisting of the iterative process based on spin selection and the iterative process based on the selection of the trotter layer, is reversed between the two.

[0159] In this embodiment, after steps 300, 302, and 304 are performed, step 700 repeats step 702 a predetermined number of times.

[0160] Step 702 also has a repeating structure. That is, step 702 includes step 720, which repeats step 722 a predetermined number of times, and step 334, which updates the temperature and magnetic field values ​​according to a schedule after the completion of step 720.

[0161] Step 722 includes step 730, which repeats the first iteration process 732 as many times as there are Trotter layers.

[0162] The first iteration 732, similar to the first iteration 622 shown in Figure 17, includes a step 350 of randomly selecting a plurality of spins and a step 752 of performing a second iteration 754 for each selected spin.

[0163] The second iteration process 754 includes the same steps 370 and 530 as shown in Figure 17, and step 372 which performs a state update process 374 for each of the spins included in the spin group determined by the processes up to step 530. In other words, in the third embodiment, the spins are selected first and then the state update process 374 is performed for each trotter layer, whereas in the fourth embodiment, the trotter layer is selected first and then the state update process 374 is performed for each spin. Therefore, in this fourth embodiment, the nesting structure is reversed compared to the third embodiment.

[0164] 5.2 Operation The processing flow in this fourth embodiment is almost the same as that shown in the third embodiment, except that the nesting structure is reversed, as described above. In the third embodiment, the process of randomly selecting spins comes first, followed by the processing of the trotter layer, whereas in the fourth embodiment, the process of selecting the trotter layer comes first, followed by the processing of randomly selecting spins and processing each spin, which is nested thereafter. It should be noted that, therefore, a detailed explanation of the operation flow of the method according to the fourth embodiment will not be repeated here. 6 Fifth Embodiment 6.1 Configuration The quantum inspired optimization methods according to the first to fourth embodiments described above relate to SA. However, other methods exist besides SA as quantum inspired optimization methods that use the Ising model. One such method is the simulated bifurcation (SB) calculation, which is introduced in Non-Patent Literature 4 and Non-Patent Literature 5. According to Non-Patent Literature 5, the SB calculation is based on a classical computer simulation of a dynamical system that exhibits bifurcation phenomena.

[0165] In SB calculations, a single oscillator is represented by a variable x representing position and a variable y representing momentum. To account for a system containing multiple oscillators, the position and momentum variables for the i-th oscillator are generally represented as xi and yi, respectively. In SB calculations, the interacting oscillators gradually diverge into either positive or negative directions over time. Simultaneously, the energy of the entire system decreases. The final variable xi obtained after a predetermined time has elapsed is considered to correspond to the i-th spin in the Ising model when binarized. This quantum-inspired optimization method can also be applied to these SB calculations.

[0166] Figure 19 shows, in flowchart form, the control structure of a computer program when this disclosure is applied to an SB calculation method, which is one embodiment of a quantum-inspired optimization method. Referring to Figure 19, this program includes step 800, which repeats the iterative process 802 a predetermined number of times.

[0167] The iterative process 802 includes a step 820 in which the spins included in the spin array are divided into multiple spin groups that can be processed in parallel (but do not have a causal relationship with each other) based on the causal relationships between the spins. In this example, the spins correspond to the position variables xi of the oscillator. This division of spin groups is performed by determining the spin groups (groups of the position variables xi of the oscillator) for a given spin in the first layer (and second layer), and then selecting spin groups from among them that do not have overlapping parts.

[0168] The program further includes step 822, which performs step 824, for each of the spin groups obtained in step 820, the calculation shown by the following equation (equation 10). The equation shown in step 824 is based on the one described in Non-Patent Document 5.

[0169] In this equation, Δyi represents the increment in momentum yi over time, Δt is the time step expressed as an integer, γ0 is a constant, Ji is the vector of the i-th row of the coupling coefficient matrix J between spins, Jij is its j-th element, and N is the number of spins.

[0170] In step 820, there are no overlapping spins (spins that influence each other) among the divided spin groups. Therefore, step 824 can be performed in parallel for each spin group.

[0171] The iteration 802 further includes step 826, after step 822 is completed, performing step 828 of the time evolution for each spin.

[0172] Step 826 includes step 840, which adds the increment Δyi calculated in step 824 to the momentum variable yi corresponding to the spin being processed (let's call it the position variable xi), and step 842, which repeats step 844 a predetermined number of times for the position variable xi and the momentum variable yi of the oscillator to be processed.

[0173] In step 844, the momentum variable yi and the position variable xi of the oscillator being processed are updated by the following equation (Equation 10).

[0174] The FX and FP appearing in this equation are expressed by the following equation (Equation 11).

[0175] Here, δt = Δt / M (where M is the substep number) is the sub-timestep, α is a control parameter that increases as a function of time t, α0 and β0 are constant parameters, η is a parameter that changes depending on α(t), and hi represents the strength of the external magnetic field.

[0176] When step 842 is completed, the processing in step 828 for the spin being processed is finished. If there are still spins to be processed, the processing in step 820 is performed, a group capable of parallel computing is selected, and the processing in step 828 is executed for those spins. Once the execution of step 828 is completed for all spins, the execution of step 826 is completed.

[0177] The iteration process 802 further includes a step 830 which is performed after step 826 and adds Δα to α to complete one execution of the iteration process 802.

[0178] Once the iterative process 802 has been executed a predetermined number of times, the spin bifurcation calculation process is terminated.

[0179] In spin bifurcation calculations, the energy value of the system decreases as time evolution progresses according to step 826, and the system is closest to its ground state when the time evolution process is complete.

[0180] According to this fifth embodiment, compared to the methods disclosed in Non-Patent Documents 4 and 5, the calculation process in step 824 in step 822 can be processed in parallel for each spin group. Therefore, as with the first to fourth embodiments, it is possible to speed up processing and improve convergence.

[0181] 6.2 Effects of the above embodiment: The effects of the first to fourth embodiments will be described in particular with reference to Figures 20 to 25.

[0182] Figure 20 is a graph showing the degree of speed improvement of the quantum inspired optimization method disclosed herein, verified by solving the MaxCut problem (G34). Figure 20 shows the results of performing the quantum inspired optimization method according to the above embodiment under various conditions, with the logarithm of the MCS on the horizontal axis and the best cut number obtained in that MCS on the vertical axis.

[0183] In Figure 20, curve 930 represents an example where only the Trotter layer is parallelized. Curve 932 represents an example where only the spins up to the first layer are used and the nesting structure is not adopted. This corresponds to the first embodiment. Curve 934 represents an example where only the spins up to the first layer are used and the nesting structure is adopted. This corresponds to the case where the third embodiment is applied. Curve 936 represents an example where the spins up to the second layer are used and the nesting structure is also adopted. This corresponds to the third embodiment. Similar results can be obtained with the fourth embodiment.

[0184] Figure 21 shows, in tabular form, the results of calculating the degree of speedup by comparing the MCS when the number of cuts reaches a predetermined ratio of BKC for each of the four embodiments shown in Figure 20. Figure 22 plots the degree of speedup of the quantum inspired optimization method according to each embodiment when the nesting structure and the number of spin layers used are changed according to the above embodiment. The degree of speedup is based on the speed when only the Trotter layer is parallelized.

[0185] In Figure 22, the black circular plots indicate the degree of speed improvement when using spins up to the first layer without employing nesting. The white triangular plots indicate the degree of speed improvement when employing a nesting structure that utilizes spins up to the first layer and performs spin selection later (i.e., using spins up to the first layer in the third embodiment). The white circular plots indicate the degree of speed improvement when employing a nesting structure that utilizes spins up to the first layer and performs spin selection first (i.e., using spins up to the first layer in the fourth embodiment). The white square plots indicate the degree of speed improvement when using spins up to the second layer and employing a nesting structure (i.e., the third embodiment).

[0186] In these calculation results, when nesting is present, the number of Trotter layers is set to 10, and partial optimization due to nesting is performed 20 times. In the comparative example (where the nested structure does not utilize the spins of the first or second layer), the number of Trotter layers is set to 200. Therefore, in the embodiments of this disclosure shown in Figures 20 to 22, the number of loops included in 1 MCS is the same as in the comparative example, except for the number of spins of the layers taken into consideration.

[0187] Referring to Figures 21 and 22, it can be seen that in the case without nesting, while there is an effect on speeding up, the degree of speed increase due to this effect is limited to about the number of spins in the first layer. Although not shown in the figures, in the processing of G34 in Gset, the number of spins in the first layer is about 5 on average, and in Figure 21, the degree of speed increase is limited to about 5.

[0188] On the other hand, when a nested structure is adopted, it is thought that partial optimization can be performed due to that nested structure. In this case, compared to curve 932 without a nested structure, it can be seen that the processing is even faster. On the other hand, these figures also show that there is no significant difference in speedup depending on whether parallelization by spin group selection is performed first or parallelization by Trotter layers is performed first. In other words, it can be predicted that the order of the two processes can be swapped to some extent, and therefore parallel processing can be predicted. For this reason, by adopting such a nested structure, it is possible to achieve a speedup equivalent to (number of Trotter layers) × (degree of parallelism of selected spin groups).

[0189] Figure 23 is a graph that shows the degree of speed improvement of the optimization method disclosed in this disclosure, when verified by solving another Maxcut problem (G9), along with the results of an optimization method that does not use causal relationships between spins.

[0190] Referring to Figure 23, curve 950 is a graph plotting the best cut number for the MCS. Curve 950 is the result when the number of Trotter layers is 150 and only one selected spin is used. Curve 950 is an example where only the Trotter layers are parallelized.

[0191] Curve 952 shows the results when the first embodiment is applied. Although only the first layer of spins was used, a nesting structure was not employed.

[0192] Curve 954 shows the results when a nesting structure like that of the third embodiment is applied. In this example, the spins used are limited to the first layer. For curve 954, the number of nesting repetitions is 50, and the number of Trotter layers is 3.

[0193] As can be seen from Figure 23, using a nested structure dramatically improves convergence compared to not using one. It is also possible to use spins up to the second layer or beyond, but this requires a larger memory capacity, so whether it is effective or not should be judged on a case-by-case basis.

[0194] Figure 24 shows the calculation results for G9, with and without nesting, with nesting and prior spin selection, and with nesting and subsequent spin selection. The degree of speed improvement is calculated by comparing the MCS when the number of cuts reaches a predetermined percentage of BKC for each of these cases.

[0195] Figure 25 plots the information shown in Figure 24. The horizontal axis of Figure 25 represents the BKC achievement rate, and the vertical axis represents the degree of speed improvement due to parallelization. In Figure 25, the black circular plots represent the case where spins up to the first layer are used without nesting (corresponding to the first embodiment). The white triangular plots represent the degree of speed improvement in the example where spins up to the first layer are used and a nesting structure is adopted (third embodiment), with spin selection performed first. The circular white plots represent the degree of speed improvement in the fourth embodiment, where spins up to the first layer are used and a nesting structure is adopted, with spin selection performed later.

[0196] In this example, the method without a nested structure uses 150 Trotter layers. For the partial optimization using a nested structure, the number of Trotter layers is set to 50 and the number of iterations to 3. As a result, the number of loops included in one MCS is the same in all of the cases shown by curves 950, 952, and 954, except for the number of spins in the layers considered.

[0197] Although not shown in the diagram, in the example where a nested structure was adopted, BKC was reached 4 times out of 100 trials.

[0198] Section 7 Consideration of the Third to Fifth Embodiments As described above, in the third and fourth embodiments, the nesting structure is reversed between the iterative processing by spin selection and the iterative processing by trotter layer selection. In the fifth embodiment, this reversal is possible to some extent, and therefore it is predicted that these parallel processing will be possible.

[0199] However, in fact, in the third and fourth embodiments, the calculation results from both are almost identical near the convergence point of the optimization calculation. This can be confirmed by performing further calculations. That is, in the nested structure including the iterative processing by spin selection and the iterative processing by Trotter layer selection described above, the timing of spin selection does not affect the result. This means that the order of the iterative processing by spin selection and the iterative processing by Trotter layer selection can be interchanged. Since the order of the two does not matter, this fact provides the basis for executing these iterative processes in parallel. Therefore, using a nested structure as in the third and fourth embodiments has the effect of enabling high-speed optimization by using hardware resources capable of parallel processing.

[0200] Further investigation into this parallel processing reveals that, for example, if only one Trotter layer is used, parallel processing cannot be expected. However, even in this case, it can be seen that optimization can be performed quickly and with good convergence. As can be seen by comparing Figure 17 and Figure 18, the computational cost of iterations per MCS is the number of Trotter layers × the number of nesting iterations. Let P be the number of Trotter layers and N be the number of nesting iterations. Let P be a value of 1, 2, 5, 10, 20, and 40. If the value of P × N is constant, the values ​​of N corresponding to these values ​​of P are 200, 100, 40, 20, 10, and 5, respectively. Figure 26 shows the results of the optimization calculations performed using these combinations.

[0201] In Figure 26, curves 760, 762, 764, 766, 768, and 770 show the results of optimization calculations when P = 1, 2, 5, 10, 20, and 40, respectively. From Figure 26, it can be seen that when P = 1, the convergence rate for the same MCS is higher than for the other P values ​​of 2, 5, 10, 20, and 40. In other words, when performing optimization calculations with the number of Trotter layers = 1, a high convergence rate can be obtained with fewer MCS compared to when using other numbers of Trotter layers. In the example shown in Figure 26, a 99% convergence rate is achieved with fewer MCS than with any other combination.

[0202] This may seem strange at first glance. However, it is important to note that Figure 26 uses MCS as the horizontal axis. This does not take into account the effect of parallel computing. Generally, when performing parallel computing, the speed can be increased by a value obtained by dividing the same amount of time by the number of parallelizations, compared to non-parallel computing. When the number of Trotter layers P is 1, the degree of parallelization = 1. If the value of P is 2, 5, 10, 20, or 40, the degree of parallelization is 2, 5, 10, 20, or 40, respectively. As a result, the computation time under each condition is reduced to MCS / degree of parallelization = MCS / P.

[0203] Taking this into consideration, Figure 27 shows the result of taking only the area to the left of the boundary 980 shown by the dotted line in Figure 26, multiplying the horizontal axis by MCS / P to consider it as the computation time during parallelization (arbitrary scale), and replotting the curves 960, 962, 964, 966, 968, and 970 from Figure 26.

[0204] In Figure 27, curve 960 shows that convergence is slower when the number of Trotter layers is 1 compared to other cases. Furthermore, it is shown that convergence slows down as the number of Trotter layers decreases. This indicates that there is little room for speedup because the contribution of parallelizing the Trotter layers cannot be expected.

[0205] Nevertheless, the fact that optimization calculations can be performed quickly and with good recording capabilities even when the number of Trotter layers is 1 is extremely significant in practical terms. This is because approximate results can be obtained in a short time, even when available computing resources are limited.

[0206] Figure 8 shows the appearance of a computer running a computer program that implements the quantum-inspired optimization method according to the above embodiment, and Figure 29 shows an example of its hardware configuration. Figures 28 and 29 are merely examples, and it goes without saying that the quantum-inspired optimization method can be implemented using computers other than those shown in Figures 28 and 29.

[0207] Referring to Figure 28, this computer system 1000 includes a computer 1050 and a keyboard 1054, a mouse 1056, and a monitor 1052, all connected to the computer 1050, for user interaction. Of course, these are just examples of configurations for when user interaction is required, and any general hardware and software available for user interaction (e.g., touch panels, voice input, pointing devices in general) can be used.

[0208] Referring to Figure 29, the computer 1050 includes a bus 1110 and, all connected to the bus 1110, a CPU (Central Processing Unit) 1100, a GPU 1102, a ROM (Read-Only Memory) 1104 that stores the computer 1050's boot-up program, a RAM (Random Access Memory) 1106 that stores instructions, system programs, and work data that make up a program when the program is executed, and a non-volatile memory SSD (Solid State Drive) 1108 connected to the bus 1110. The SSD 1108 is for storing programs executed by the CPU 1100 and GPU 1102, as well as data used by programs executed by the CPU 1100 and GPU 1102. These programs include programs for implementing the quantum-inspired optimization method described above.

[0209] The computer 1050 further includes a network interface 1118 that provides connectivity to the network 1060, and a USB port 1116 that allows for the insertion and removal of a USB memory 1058 and provides communication between the USB memory 1058 and various parts within the computer 1050. In this embodiment, the network 1060 is the internet.

[0210] The computer 1050 further includes an input / output interface 1114 connected to the bus 1110 for inputting input from the keyboard 1054 and mouse 1056 to the CPU 1100, and a display control unit 1112 connected to the bus 1110 for controlling the display on the monitor 1052 according to the control from the CPU 1100.

[0211] The computer program for operating this computer system 1000 to perform the functions described later is stored in the USB memory 1058. The USB memory 1058 is plugged into the USB port 1116, and the program is transferred to the SSD 1108. Alternatively, the program may be transmitted to the computer 1050 via the network 1060 and stored in the SSD 1108. The program is loaded into the RAM 1106 when it is executed. Part of the program may be stored in the ROM 1104.

[0212] The CPU 1100 reads the program from RAM 1106 according to the address indicated by an internal register called the program counter (not shown), interprets the instructions, reads the data necessary for executing the instructions from RAM 1106, SSD 1108, or other devices according to the address specified by the instructions, and executes the processing specified by the instructions. The CPU 1100 stores the execution result data at an address specified by the program, such as RAM 1106, SSD 1108, or a register within the CPU 1100. At this time, the value of the program counter is also updated by the program. The computer program may be loaded directly into RAM 1106 from USB memory 1058 or via a network. Among the programs executed by the CPU 1100, tasks that can be executed in parallel among the iterative processes (calculations for each Trotter layer, state updates for each selected spin, and calculations for both) can be executed by the GPU 1102 according to the instructions included in the program or according to the analysis results when the CPU 1100 executes the instructions. This makes it possible to output the results of the combinatorial optimization problem at high speed according to the quantum-inspired optimization method described above.

[0213] The program that, in cooperation with the computer 1050, implements each step of the quantum-inspired optimization method, each function of the quantum-inspired optimization apparatus, and the method for solving the path optimization problem, includes a plurality of instructions that are written and arranged to operate the computer 1050 to implement those functions. The operation method of the computer 1050 that executes these instructions is well known. Therefore, the explanation of that operation method will not be repeated here.

[0214] As described above, each embodiment of this disclosure enables the implementation of a quantum-inspired optimization method that can reduce computation time compared to conventional methods. Furthermore, this disclosure provides a fast approximate solution method for problems that require enormous computational resources, such as path optimization problems.

[0215] Each process (each function) of each embodiment described above is performed by a processing circuit (Circuitry) including one or more processors. The processing circuit may consist of one or more memories, various analog circuits, various digital circuits, and other integrated circuits in addition to the one or more processors. The one or more memories store programs (instructions) that cause the one or more processors to execute each of the processes. The one or more processors may execute each of the processes according to the programs read from the one or more memories, or they may execute each of the processes according to logic circuits that have been pre-designed to execute each of the processes. The above-mentioned processor may be various processors suitable for computer control, such as a CPU, GPU, DSP (Digital Signal Processor), FPGA (Field-Programmable Gate Array), ASIC (Application Specific Integrated Circuit), and DPU (Data Processing Unit). Multiple physically separated processors may cooperate with each other to perform the above-mentioned processes. For example, the above-mentioned processors installed in each of multiple physically separated computers may cooperate with each other via a network such as a LAN (Local Area Network), WAN (Wide Area Network), or the Internet to perform the above-mentioned processes. The above program may be installed in the memory via the network from an external server device, or it may be distributed on a recording medium such as a CD-ROM (Compact Disc Read-Only Memory), DVD-ROM (Digital Versatile Disc Read-Only Memory), or semiconductor memory, and then installed in the memory from the recording medium.

[0216] The embodiments described above relate to a two-dimensional Ising model. However, this disclosure is not limited to such embodiments. It can also be applied to a three-dimensional Ising model, and in principle, to any calculation method that uses the Suzuki-Trotter decomposition to perform Monte Carlo calculations.

[0217] With respect to the embodiments described above, the following additional notes are disclosed: (1) An optimization method according to one aspect of this disclosure is an optimization method for searching for a minimum energy state that gives the minimum energy of a spin system including a plurality of spins, the method comprising the steps of: a computer setting an initial value for the minimum energy; and a computer searching for the minimum energy state of the spin system when the spin state of the spins included in the spin system is changed, wherein the computer may search for the minimum energy state by dividing the plurality of spins into a plurality of spin groups based on the causal relationship between the plurality of spins, calculating the energy of the spin system by changing the spin state of the spins included in each spin group, and updating the minimum energy in response to obtaining an energy smaller than the value of the minimum energy.

[0218] (2) In (1) above, the search step may include a first iteration step in which the computer repeats the first process a first predetermined number of times, and a step in which the computer identifies the state of the spin system that gives the minimum energy when the first iteration step is completed as the minimum energy state.

[0219] (3) In (2) above, the first process may include a second iteration step in which the computer repeats the second process as many times as there are Trotter layers in the spin system, and a step in which the computer updates the environment of the spin system according to a predetermined update method after the completion of the second iteration step.

[0220] (4) In (3) above, the second predetermined number of times may be the number of trotter layers in the spin system.

[0221] (5) In (3) or (4) above, the second process may include the steps of: the computer dividing the spins included in the spin system into a plurality of spin groups based on the causal relationships between the spins included in the spin system; and the computer updating the minimum energy value for each spin included in the plurality of spin groups in response to the fact that the energy value of the spin system when the spin state of the spin is changed has become smaller than the minimum energy value.

[0222] (6) In (5) above, the step of dividing the spins included in the spin system may include the step of identifying a spin group that includes the spin and the spins affected by the change in the spin state of the spins in the spin system when the computer changes the spin state of the spins.

[0223] (7) In (6) above, the step of identifying the spin group may include a selection step in which the computer selects one or more arbitrary spins of the spin system as target spins; a step in which, for each of the target spins selected in the selection step, the computer determines which spins in the spin system are directly affected by the target spin; and a step in which the computer identifies the group including the target spins and the spins determined in the determination step as the spin group.

[0224] (8) In (7) above, the determination step may include a determination step in which the computer determines the spins that are directly affected by the target spin from the matrix elements of the Hamiltonian matrix describing the spin system that have non-zero values.

[0225] (9) In (8) above, the determination step may further include the step in which the computer, in the spin system, takes the spin determined in the determination step as the target spin and determines the spins that are directly affected by the target spin from the matrix elements of the Hamiltonian matrix that have non-zero values.

[0226] (10) In (7) above, the determination step may further include the step in which the computer, in the spin system, takes the spin determined in the determination step as the target spin and determines the spins that are directly affected by the target spin from the matrix elements having non-zero values ​​of the Hamiltonian matrix describing the spin system.

[0227] (11) In (3) above, the spin system may include one or more Trotter layers, and the second process may include the steps of the computer selecting spins based on the causal relationships between the spins included in the spin form and running a loop until all spins are selected, and the third iteration step in which the computer repeats the third process for each of the one or more Trotter layers as many times as there are Trotter layers.

[0228] (12) In (3) above, the spin system may include one or more Trotter layers, the second process may include a third iteration step in which the computer repeats the third process as many times as there are Trotter layers, the third process may include a step in which the computer divides the spins included in the spin system into a plurality of spin groups, each containing a first predetermined number of spins, based on the causal relationships between the spins included in the spin system, and an update step in which the computer, for each spin included in the plurality of spin groups, updates the minimum energy value with the energy value of the spin system when the spin state of the spin is changed, in response to the energy value of the spin system when the spin state of the spin is changed becoming smaller than the minimum energy value.

[0229] (13) In (3) above, the spin system may include one or more Trotter layers, and the second process may include a fourth iteration step in which the computer repeats for each of the one or more Trotter layers, and a fifth iteration step in which the computer selects spins based on the causal relationships between the spins included in the spin form and iterates through the loop until all spins are selected.

[0230] (14) In any one of (3) to (13) above, the update step may include a selective update step in which the computer selectively performs a process of updating the minimum energy state each time the energy of the spin system when the spin state of each spin included in the spin group is changed in the one or more Trotter layers is less than the minimum energy value calculated in the one or more Trotter layers, and a process of not updating the minimum energy state.

[0231] (15) A computer program relating to one aspect of this disclosure causes a computer to execute any one of the optimization methods described in (1) to (13) above.

[0232] (16) An optimization apparatus according to one aspect of this disclosure includes a computer programmed to perform any one of the optimization methods described in (1) to (13) above.

[0233] (17) A method for solving a path optimization problem according to one aspect of this disclosure includes the steps of: a computer converting a path optimization problem formulated in a predetermined format into an energy problem of a spin system; a computer determining a solution for the spin system by applying any one of the optimization methods (1) to (13) above to the spin arrangement of the spin system; and a computer converting the solution into a solution for the path optimization problem.

[0234] The embodiments disclosed herein should be considered in all respects to be illustrative and not restrictive. The scope of this disclosure is not indicated by the description in the detailed description of the disclosure but by each claim, and all modifications within the meaning and scope equivalent to the wording of the claims are intended to be included.

[0235] 50, 90 Ising model 60, 62, 64 Nodes 70, 72, 240, 242, 244, 246, 470, 472, 474 Spins 100 Trotter layer 158 Spin update table 180 Graph 190, 192, 230, 460 Spin group 194 Edge set 220, 222 Vector 224 Hamiltonian 226 Connection matrix 374 State update process 380, 382, ​​430, 432, 580, 582, 584, 930, 932, 934, 936, 950, 952, 954, 960, 962, 964, 966, 968, 970 Curve 384 Line 445 Spin group table 622, 732 First iteration process 682, 754 Second Iteration 802 Iteration 980 Boundary 1000 Computer System 1050 Computer 1052 Monitor 1054 Keyboard 1056 Mouse 1058 USB Memory 1060 Network 1100 CPU 1102 GPU 1104 ROM 1106 RAM 1108 SSD 1110 Bus 1112 Display Control Unit 1114 Input / Output I / F 1116 USB Port

Claims

1. An optimization method for searching for a minimum energy state that gives the minimum energy of a spin system, comprising the step of a computer searching for the minimum energy of the spin system when the spin state of the spins included in the spin system is changed, wherein the searching step is an optimization method in which the computer divides the spins of the spin system into a plurality of spin groups based on the causal relationships between the spins, and for each of the spin groups, searches for the minimum energy by calculating the energy of the spin system by changing the spin state of the spins included in the spin group.

2. The optimization method according to claim 1, wherein the search step includes: a step of identifying a spin group including a spin and a spin affected by the change in the spin state of the spin when the spin state of the spin system is changed by the computer; a state update step in which the computer selectively performs, for each of the spin groups, a process of updating the minimum energy state with respect to the spins included in the spin group, according to whether the energy of the spin system when the spin state of the spin is changed is less than the minimum energy value, and a process of not updating the minimum energy state; an environment update step in which the computer updates the environment of the spin system in response to the completion of the state update step; and an iterative execution step in which the computer repeatedly performs the first process including the step of identifying the spin group, the state update step, and the environment update step until a predetermined termination condition is met.

3. The optimization method according to claim 2, wherein the step of identifying the spin group comprises: a selection step in which a computer selects one or more arbitrary spins of the spin system as target spins; a step in which, for each of the target spins selected in the selection step, the computer determines the spins in the spin system that are directly affected by the target spin; and a step in which the computer identifies the group including the target spins and the spins determined in the determination step as the spin group.

4. The optimization method according to claim 2 or 3, wherein the determining step includes a determination step in which a computer determines the spins directly affected by the target spin from the matrix elements having non-zero values ​​of the Hamiltonian matrix describing the spin system.

5. The optimization method according to claim 4, wherein the determination step further includes the step of a computer determining, in the spin system, the spin determined in the determination step to be the target spin, and determining the spins directly affected by the target spin from the matrix elements of the Hamiltonian matrix that have non-zero values.

6. The optimization method according to any one of claims 2 to 5, wherein the first process includes a step in which a computer performs a predetermined second process a predetermined number of times, the second process includes a first iteration step in which the computer repeatedly performs the step of identifying the spin group and the state update step and terminates in response to a predetermined first termination condition being met, and a second iteration step in which the computer repeatedly performs the first iteration step and the environment update step and terminates in response to a termination condition being met.

7. The optimization method according to claim 6, wherein the first iteration step includes a step in which a computer repeatedly performs the step of identifying the spin group and the state update step, and terminates in response to the number of iterations reaching a predetermined number.

8. The optimization method according to claim 6 or 7, wherein the second iteration step includes a step in which the computer repeatedly performs the first iteration step and the environment update step, and terminates in response to the number of iterations of the first iteration step and the environment update step reaching a predetermined number of iterations.

9. The optimization method according to any one of claims 2 to 8, wherein the spin system comprises a plurality of trotter layers, and the state update step includes a trotter layer-specific state update step in which a computer selectively performs a process of updating the minimum energy state and a process of not updating the minimum energy state, depending on whether the energy of the spin system when the spin state of each spin included in the spin group is changed in the plurality of trotter layers is less than the minimum energy value calculated by the plurality of trotter layers.

10. The optimization method according to any one of claims 2 to 8, wherein the spin system comprises a single Trotter layer, and the state update step includes a step of selectively performing a process of updating the minimum energy state and a process of not updating the minimum energy state, depending on whether the energy of the spin system when the spin state of each spin included in the spin group is changed in the single Trotter layer is less than the minimum energy value calculated by the single Trotter layer.

11. The optimization method according to claim 6, wherein the first iteration step includes, in each Trotter layer of the spin system, the step of a computer identifying the spin group, and the step of a computer performing the state update step for each spin included in the spin group in all Trotter layers in the step of identifying the spin group.

12. The optimization method according to any one of claims 2 to 5, wherein the spin system comprises one or more Trotter layers, the first process includes a second iteration step in which a computer performs a predetermined second process a predetermined number of times, the second process includes a third iteration step in which a computer performs a predetermined third process a number of times corresponding to the number of Trotter layers, and a step in which the computer performs the environment update step after the completion of the third iteration step, the third process includes a fourth iteration step in which the computer repeatedly performs the steps of identifying the spin group and the state update step and terminates in response to a predetermined first termination condition being met, and a fifth iteration step in which the computer repeatedly performs the fourth iteration step and the environment update step and terminates in response to a termination condition being met.

13. The optimization method according to claim 12, wherein the state update step includes a step in which the computer selectively performs a process of updating the minimum energy state and a process of not updating the minimum energy state, depending on whether the energy of the spin system when the spin state of each spin included in the spin group is changed in the one or more Trotter layers is less than the minimum energy value calculated by the one or more Trotter layers.

14. A computer program that causes a computer to execute the optimization method described in any one of claims 1 to 13.

15. An optimization apparatus including a computer programmed to perform the optimization method according to any one of claims 1 to 13.

16. A method for solving a path optimization problem, comprising the steps of: a computer converting a path optimization problem formulated in a predetermined format into an energy problem of a spin system; a computer determining a solution for the spin system by applying the optimization method described in any one of claims 1 to 13 to the spin arrangement of the spin system; and a computer converting the solution into a solution for the path optimization problem.

Images