Monte Carlo path integral-simulation quantum annealing method and system based on FPGA (Field Programmable Gate Array)

Abstract

The invention discloses a Monte Carlo path integral-simulated quantum annealing method realized based on an FPGA (Field Programmable Gate Array), which mainly solves the problem of low running time efficiency of the existing Monte Carlo path integral-simulated quantum annealing algorithm on a CPU (Central Processing Unit), and the implementation scheme comprises the following steps: initializing an original matrix in a Bram space; the FPGA generates a random number to perform random operation on the original matrix to obtain a new matrix; calculating potential energy item difference values, kinetic energy item difference values and Hamiltonian difference values of the new matrix and the original matrix in parallel; judging whether to update the original matrix; selecting a matrix with the minimum potential energy item in the original matrix after the updating operation is executed as a current optimal scheme solution; and updating the coupling coefficient, and continuously carrying out iteration until a set number of iterations, so as to obtain an optimal solution, namely a final optimal path. According to the method, the operation time efficiency is greatly improved, the operation time is improved by 79 times compared with that of a traditional Monte Carlo path integral-simulated quantum annealing algorithm, and the method can be used for unmanned aerial vehicle routing planning in a scene with a high real-time requirement.

Description

Technical Field

[0001] The present invention belongs to the field of physical technologies and relates to a Monte Carlo path integral-simulated quantum annealing method and system, which can be used for combinatorial optimization such as unmanned aerial vehicle routing planning and express delivery schemes in scenarios with high real-time requirements. Background Art

[0002] Quantum mechanics has been regarded as one of the most revolutionary developments in the history of science for more than 100 years. Quantum computing can achieve exponentially high-speed operations by utilizing quantum mechanical concepts, bringing a major breakthrough in the field of computing and becoming one of the important directions for future scientific development. With the continuous in-depth exploration of the field of quantum computing, people have tried to combine quantum computing with intelligent optimization algorithms to form quantum optimization algorithms. The quantum annealing algorithm QA is an algorithm commonly used to solve combinatorial optimization problems. It can find a relatively optimal solution of the system within a specified time. By utilizing the quantum tunneling effect, it can avoid being trapped in the dilemma of local minima. And many difficult problems in real life can be transformed into decision-making problems and then attributed to the task of combinatorial optimization solution, such as linear programming, path optimization, minimum and maximum problems, etc. With the progress of quantum computing technology, the quantum annealing algorithm has attracted much attention and research. In practical applications, it has made significant breakthroughs. For example, it has planned the operation routes of 10,000 taxis for Volkswagen to solve the traffic flow optimization problem. In addition, the quantum annealing algorithm has also been successfully applied to many fields such as training deep neural networks and image-based fault detection and diagnosis.

[0003] By means of Monte Carlo path integral, a D-dimensional quantum system can be mapped into a (D + 1)-dimensional classical physical model, enabling the performance of the quantum annealing algorithm to be simulated on a classical computer using P copies of the classical Hamiltonian. However, even if the algorithm is implemented on a classical computer, it is difficult to meet the requirements of low latency and real-time performance. Some scholars have begun to search for methods to accelerate the operation of the PIMC-QA algorithm using hardware, including ASIC, GPU, and FPGA. But there are corresponding problems. For example, the reconfigurability of ASIC is relatively poor, and the GPU is restricted by a fixed hardware structure, with relatively low flexibility and unable to maximize the algorithm efficiency. At the same time, the price cost is relatively high. Indian scholars used FPGA to accelerate the Ising model, but it does not have generalization for combinatorial optimization problems in real life. Japanese scholars used FPGA to accelerate the PIMC-QA algorithm to solve common combinatorial optimization problems, but there are problems such as large consumption of lookup tables, flip-flops, and Bram storage space, and there is still room for improvement in the acceleration effect. Summary of the Invention

[0004] The object of the present invention is to propose a method and system for implementing Monte Carlo path integral - simulated quantum annealing based on FPGA in view of the deficiencies of the above - mentioned existing technologies, so as to further improve the operation efficiency of Monte Carlo path integral - simulated quantum annealing and reduce the resource overhead on FPGA by designing a general - purpose hardware architecture for the PIMC - QA algorithm.

[0005] To achieve the above object, the technical solution of the present invention includes:

[0006] Technical solution 1: The Monte Carlo path integral - simulated quantum annealing method implemented based on FPGA provided by the present invention is characterized in that a quantum annealing machine is simulated on the FPGA for path planning, which includes:

[0007] (1) Initialize the number of replicas P and the coupling coefficient J in the Monte Carlo path integral - simulated quantum annealing algorithm, set the initial number of iterations M and the position matrix d of the path integral Monte Carlo quantum annealing algorithm; g

[0008] (2) Initialize P quadratic unconstrained binary optimization matrices QUBO to form the original matrix set S, and store them in P BRAM1s in the FPGA;

[0009] (3) Read the original matrix set S from the BRAM1, generate random numbers through the FPGA, and use the random numbers to perform random operations on the original matrix set S to update and obtain P new QUBO matrices, forming a new matrix set S′, and store them in another P BRAM2s in the FPGA;

[0010] (4) Calculate the potential energy term difference ΔH p , kinetic energy term difference ΔH k , and classical Hamiltonian difference ΔH c between the original matrix set S and the new matrix set S′ respectively;

[0011] (5) Compare the potential energy term difference ΔH p and the classical Hamiltonian difference ΔH c . If one of the differences is less than 0, or the classical Hamiltonian difference ΔH c satisfies the Monte Carlo criterion, then replace the corresponding matrix in the original matrix set S with the corresponding matrix in the new matrix set S′, and select the matrix with the minimum potential energy term in the updated original matrix set S to obtain the optimal solution S of the current Monte Carlo cycle; best

[0012] (6) Update the coupling coefficient J g , subtract 1 from the number of iterations M, and judge whether the current number of iterations M is greater than 0:

[0013] If so, return to step (3);

[0014] Otherwise, the FPGA stops working and the optimal path result, i.e., the matrix S of the current iteration, is obtained. best 。

[0015] Technical solution 2: A Monte Carlo path integral - simulated quantum annealing system implemented based on an FPGA provided by the present invention includes:

[0016] A Bram storage space for storing the original matrix set and the new matrix set;

[0017] A Bram space access module for controlling the read - write addresses of the Bram storage space;

[0018] A random number emission module for generating random numbers between 0 and N;

[0019] A random operation module for performing random operations on the original matrix according to the random numbers to obtain a new matrix;

[0020] A potential energy term calculation module for solving the potential energy terms of the original matrix and the new matrix;

[0021] A kinetic energy term calculation module for solving the kinetic energy terms of the original matrix and the new matrix;

[0022] A judgment and update module for calculating the potential energy term interpolation, kinetic energy term difference, and classical Hamiltonian difference, and judging whether to replace the old matrix with the new matrix according to these results.

[0023] An iteration counting module for, after the update module finishes execution, subtracting 1 from the iteration number M and judging whether the iteration number M is greater than 0: if so, returning to the random operation module; otherwise, the FPGA stops working.

[0024] Technical solution 3: An electronic device based on an FPGA provided by the present invention includes a memory, a reconfigurable logic operation unit, and a computer program running on the FPGA, characterized in that when the program is executed on the FPGA, any one of the Monte Carlo path integral - simulated quantum annealing in Technical solution 1 and Technical solution 2 is implemented.

[0025] Compared with the prior art, the present invention has the following advantages:

[0026] First, in Technical solution 1 of the present invention, random integers that meet the Monte Carlo randomness requirements are generated by the FPGA, eliminating the need for other external devices to generate random numbers and transfer them to the FPGA through a high - speed bus, reducing the time consumed by the transfer of random numbers between memories and further improving the operating efficiency of the system.

[0027] Second, in Technical solution 1 of the present invention, the potential energy term difference ΔH is calculated on the FPGA p , the kinetic energy term difference ΔH kThe difference ΔH between the classical Hamiltonian c , for the difference ΔH of the potential energy term p and the difference ΔH of the classical Hamiltonian c are compared to determine whether to replace the corresponding matrix in the original matrix set S with the corresponding matrix in the new matrix set S'. The above operation can prevent the optimal path result S best from falling into a local minimum.

[0028] Thirdly, all modules in Technical Solution 2 of the present invention can operate at a 200M clock. At the same time, the output of the random operation module is connected to the inputs of the potential energy term calculation module and the kinetic energy term calculation module to calculate the kinetic energy term and the potential energy term in parallel, giving full play to the advantage of parallel execution. Compared with other implementation methods, the present invention has a faster operation rate. BRIEF DESCRIPTION OF THE DRAWINGS

[0029] Figure 1 is the implementation flowchart of Embodiment 1 of the present invention;

[0030] Figure 2 is Figure 1 the schematic flowchart of generating random numbers by FPGA in

[0031] Figure 3 is the structural block diagram of Embodiment 2 of the present invention;

[0032] Figure 4 is the optimization effect diagram of the present invention;

[0033] Figure 5 is the FPGA time consumption diagram of the present invention. DETAILED DESCRIPTION OF THE EMBODIMENTS

[0034] In order to enable those skilled in the art to better understand the solution of the present invention, the following will clearly and completely describe the technical solution in the embodiments of the present invention with reference to the accompanying drawings in the embodiments of the present invention. Obviously, the described embodiments are only a part of the embodiments of the present invention, rather than all the embodiments. Based on the embodiments of the present invention, all other embodiments obtained by those of ordinary skill in the art without creative efforts shall fall within the protection scope of the present invention.

[0035] Embodiment 1, a Monte Carlo path integral-simulated quantum annealing method based on FPGA.

[0036] Referring to Figure 1 . This example is implemented based on FPGA, that is, a quantum annealing machine is simulated on the FPGA for path planning. The implementation steps are as follows:

[0037] Step 1, setting of initial parameters:

[0038] Set the number of replicas P, temperature T, magnetic field strength Γ, initial number of iterations M, maximum number of iterations Max, and position matrix d in the Carlo path integral - simulated quantum annealing algorithm.

[0039] In this example, it is set but not limited to: P = 8, T = 50, Max = 10000, M = 10000, Γ = 480.

[0040] Step 2, Initialize the original matrix set E and calculate the coupling coefficient J g .

[0041] (2.1) Set a matrix of all 0s with N rows and N columns, where the rows represent the node numbers in the path planning and the columns represent the time points;

[0042] (2.2) In the current row, randomly select an integer a between 0 and N, and assign the element in the a-th column of the current row to 1;

[0043] (2.3) Repeat step (2.2) until all rows are selected to obtain an original matrix s;

[0044] (2.4) Repeat (2.1) to (2.3) until all P matrices are initialized to form the original matrix set E,

[0045] and store it in P Brams on the FPGA;

[0046] E = {s 1 , s 2 , … s m …, s P},

[0047] where s m is the m-th matrix in the original matrix set E, and m is an integer belonging to 1 to P; each matrix is a randomly arranged 0, 1 matrix;

[0048] (2.5) Calculate the coupling coefficient J g :

[0049] J g = sin(-T / 2) * log(tanh(Γ / PT)).

[0050] where T is the temperature, Γ is the magnetic field strength, and P is the number of replicas;

[0051] Step 3, The FPGA generates random numbers.

[0052] Refer to Figure 2 , The implementation of this step includes the following:

[0053] (3.1) Design an oscillator including three cascaded NOT gates on the FPGA backplane;

[0054] (3.2) Connect two oscillators pairwise to form a set of oscillation loops. At the moment of the sampling clock, connect the output of the first oscillator a to the input of the second oscillator b, and simultaneously connect the output of the second oscillator b to the input of the first oscillator a;

[0055] (3.3) Arrange 19 oscillation loops in a parallel up-and-down manner in the FPGA, perform an exclusive OR operation on all the output results, and output a single-bit data;

[0056] (3.4) Use 32 registers in the FPGA to collect the single-bit data for 32 clock cycles to obtain the initial random number seed, which serves as the initial state of the m-sequence;

[0057] (3.5) According to the primitive polynomial of the m-sequence: x 32 +x 23 +x 22 +x 22 +1, predict the m-sequence value after log 2 (N) clock cycles;

[0058] (3.6) Intercept the lower log 2 (N) bits of the newly generated m-sequence to obtain a random number that meets the requirements, that is, generate a random integer between 0 - N that satisfies the Monte Carlo randomness requirement in the FPGA.

[0059] Step 4, perform random operations on the original matrix S to obtain a new matrix set E'.

[0060] (4.1) Take out two different random numbers c and f from the random numbers;

[0061] (4.2) Take out a matrix from the original matrix set E, and judge whether the second random number f is odd.

[0062] If so, the current matrix executes (4.3).

[0063] Otherwise, execute (4.4);

[0064] (4.3) Perform a random operation to symmetrically flip all the columns between c and f;

[0065] (4.4) Perform a random operation, use c and f as row indices, and swap the positions of the two rows;

[0066] (4.5) Repeat steps (4.1) to (4.4) until all P matrices have completed random operations to obtain a new matrix set E',

[0067] E' = {s ′1 , s ′2 , … s ′m …, s′P},

[0068] where s ′m represents the m-th new matrix s' in the new matrix set E', and m is an integer belonging to 1 to P; each matrix is a 0, 1 matrix randomly arranged with N rows and N columns.

[0069] Step 5, calculate the potential energy terms of the original matrix and the new matrix and their differences.

[0070] (5.1) Calculate the potential energy term of the original matrix

[0071]

[0072] where d i,j is the element at the abscissa i and ordinate j in the position matrix d, representing the distance between node i and node j;

[0073] represents the element at the abscissa i and ordinate k in s m ;

[0074] s mT is the transition matrix obtained by shifting s m one bit to the left, represents the element at the abscissa k and ordinate j in the transition matrix s mT ;

[0075] represents the element at the abscissa i and ordinate j in s m ;

[0076] W max represents the element with the largest value in the position matrix d;

[0077] (5.2) Calculate the potential energy term of the new matrix

[0078]

[0079] where s ′m is the m-th new matrix in the new matrix set E' obtained after random operation;

[0080] represents the element at the abscissa i and ordinate k in s ′m ;

[0081] s ′mT is the new transition matrix obtained by shifting s ′m one bit to the left;

[0082] is s ′mTThe element with abscissa k and ordinate j;

[0083] Denote s ′m The element with abscissa i and ordinate j in

[0084] W max Denote the element with the largest value in the position matrix d.

[0085] (5.3) According to and Obtain the potential energy term difference ΔH p :

[0086]

[0087] Step 6, calculate the kinetic energy terms and their differences of the original matrix and the new matrix.

[0088] (6.1) Calculate the kinetic energy term of the original matrix

[0089]

[0090] Among them, Denote the element with abscissa i and ordinate j in the (m - 2)-th matrix in the original matrix set E;

[0091] Denote the element with abscissa i and ordinate j in the m-th matrix in the original matrix set E;

[0092] Denote the element with abscissa i and ordinate j in the (m + 2)-th matrix in the original matrix set E;

[0093] (6.2) Calculate the kinetic energy term of the new matrix

[0094]

[0095] Among them, Denote the element with abscissa i and ordinate j in the (m - 1)-th new matrix s ′m-1 in the new matrix set; Denote the element with abscissa i and ordinate j in the m-th new matrix s ′m in the new matrix set;

[0096] Denote the element with abscissa i and ordinate j in the (m + 1)-th new matrix s ′m+1 in the new matrix set;

[0097] (6.3) According to and Calculate the kinetic energy term difference ΔH k :

[0098]

[0099] Step 7, calculate the classical Hamiltonian difference ΔH c .

[0100] According to the kinetic energy term difference ΔH k and the potential energy term difference ΔH p calculate the classical Hamiltonian difference ΔH c :

[0101] ΔH c = ΔH p / P - J g ×ΔH k .

[0102] where P is the number of replicas and J g is the coupling coefficient.

[0103] Step 8, obtain the current optimal path result through the classical Hamiltonian difference ΔH c .

[0104] (8.1) Judge whether to replace the old matrix S with the new matrix S′ according to ΔH c :

[0105] If one of the two parameters, the potential energy term difference ΔH p and the classical Hamiltonian difference ΔH c is less than 0, or the classical Hamiltonian difference ΔH c satisfies the Monte Carlo criterion, then replace the original matrix S with the new matrix S′ to obtain the currently updated matrix S″;

[0106] (8.2) Compare the potential energy terms of the updated matrix S″ with the other unupdated matrix potential energy terms in the original matrix set E, and take the minimum value as the optimal matrix S for the current iteration best .

[0107] Step 9, update the iteration number M to obtain the optimal task path.

[0108] Update the magnetic field strength: Γ = Γ - Γ / Max;

[0109] Update the iteration number: M = M - 1, and then judge whether the current iteration number M is greater than 0;

[0110] If so, return to Step 2;

[0111] Otherwise, take the optimal matrix S best after the current iteration as the optimal task path for path planning.

[0112] Embodiment 2. Monte Carlo Path Integration - Simulated Quantum Annealing System Implemented Based on FPGA

[0113] Referring to Figure 3 , this example includes: Bram storage space 1, Bram space access module 2, random number emission module 3, random operation module 4, potential energy term calculation module 5, kinetic energy term calculation module 6, judgment and update module 7, and iteration counting module 8. Its working principle is as follows:

[0114] The random number emission module 3 generates random integers between 0 and N. After the random operation module 4 takes out two different random numbers from the random number emission module 3, it sends a request signal to the Bram space access module 2. After receiving this request signal, the Bram space access module 2 reads the original matrix set data from the Bram storage space 1 and sends it to the random operation module 4;

[0115] The random operation module 4 inputs the obtained new matrix set into the Bram space access module 2, potential energy term calculation module 5, and kinetic energy term calculation module 6 simultaneously. The potential energy term calculation module 5 and kinetic energy term calculation module 6 calculate the potential energy term and kinetic energy term of the new matrix and the original matrix in parallel, and complete the calculation of the kinetic energy term and potential energy term successively;

[0116] The potential energy term calculation module 5 sends a start signal to the judgment and update module 7. The judgment and update module 7 starts to calculate the potential energy term interpolation, kinetic energy term difference, and classical Hamiltonian difference, and judges whether to replace the old matrix with the new matrix based on these results. If replacement is required, the judgment and update module 7 sends a request signal to the Bram space access module 2 to read the corresponding matrix S' in the new matrix set E' from the Bram storage space 1 and replace the corresponding matrix S in the original matrix set E.

[0117] After the judgment and update module 7 finishes the update, the iteration counting module 8 subtracts 1 from the iteration number M and judges whether the iteration number M is greater than 0: If so, it returns to the random operation module 4 to continue working; otherwise, the FPGA stops working and takes the optimal matrix S after the current iteration best As the optimal task path for path planning.

[0118] Embodiment 3. Electronic device based on FPGA.

[0119] The electronic device with FPGA provided in this example includes, but is not limited to, a memory, a reconfigurable logic operation unit, and a computer program running on the FPGA.

[0120] The memory includes common DDR, Bram, flip - flops, and configuration memory SRAM in the FPGA, and is used to store the data generated in the Monte Carlo path integration - simulated quantum annealing algorithm and the bit - stream information generated by the EDA tool;

[0121] The reconfigurable logic operation unit is the underlying circuit structure in the FPGA. The EDA development tool will convert the computer program into the underlying reconfigurable logic operation unit in the FPGA to implement the functions required by the Monte Carlo path integral - simulated quantum annealing algorithm, which includes: reading data from the memory; calculating the kinetic energy term and potential energy term of the original matrix and the new matrix in parallel; determining whether the original matrix needs to be updated; and obtaining the optimal solution matrix after iterative counting.

[0122] The computer program running on the FPGA is used to implement any one of the Monte Carlo path integral - simulated quantum annealing steps in Embodiment 1 and Embodiment 2 above.

[0123] In addition, if the bit - stream information generated by the EDA tool is stored in the Flash memory on the FPGA, this FPGA electronic device can be used as an independent device and placed separately in the edge device, and the optimal solution matrix is communicated with the edge device through other high - speed buses such as PCIE, network interface, and fiber - optic interface.

[0124] The effects of the present invention can be further illustrated by the following simulation results:

[0125] I. Simulation conditions

[0126] The software platform for the simulation experiment of the present invention is: Windows 11 operating system, Matlab 2018, Vivado2018.3.

[0127] Simulation parameters: P = 8, T = 50, Max = 10000, M = 10000, Γ = 480.

[0128] II. Simulation content

[0129] Simulation 1, under the above - mentioned simulation conditions, use the method of the present invention to simulate the quantum annealing machine on the FPGA for path planning, and use the optimal matrix S of each iteration best as the optimal task path. The potential energy term of the optimal matrix S best represents the total distance of the optimal task path. The results are as Figure 4 follows:

[0130] From Figure 4 it can be seen that the total distance of the initial optimal task path is 9300482. As the number of iterations increases, the total distance of the optimal path gradually decreases. When the number of iterations reaches MAX, the total distance of the final optimal path decreases to 896570, and the optimization effect relative to the initial optimal task path reaches more than 85%.

[0131] Simulation 2. Under the above simulation conditions, the path integral Monte Carlo-simulated quantum annealing algorithm implemented by the method of the present invention and the traditional i5-12500h CPU are compared in terms of operation time, and the results are as follows Figure 5 , where 5(a) represents the time consumed for one iteration of the path integral Monte Carlo-simulated quantum annealing algorithm on the FPGA in the present invention, and the unit is microseconds; 5(b) represents the time consumed for implementing the path integral Monte Carlo-simulated quantum annealing algorithm on the CPU.

[0132] From Figure 5 (a), it can be seen that the present invention only needs 2.3025 microseconds to obtain the result of one iteration on the FPGA, that is, the final optimal path can be obtained by running for 23.025 milliseconds.

[0133] From Figure 5 (b), it can be seen that the CPU needs 1.9 seconds to obtain the final optimal path.

[0134] Comparing Figure 5 (a) and Figure 5 (b), it can be known that the path integral Monte Carlo-simulated quantum annealing algorithm implemented by the present invention based on the FPGA shortens the running time by more than 79 times compared with the CPU, greatly improving the time efficiency of path optimization.

[0135] It should be noted that the step numbers in the specification and claims of the present invention are only for clearly describing the implementation embodiments of the present invention for easy understanding, and their serial numbers are not limited.

Claims

1. A Monte Carlo path integral-simulated quantum annealing method based on FPGA, characterized in that: exist The quantum annealing machine is simulated on FPGA for path planning, which includes: (1) Initialize the number of replicas P, temperature T, magnetic field strength Γ in the Monte Carlo path integral-simulated quantum annealing algorithm, and set the initial number of iterations M and position matrix d of the path integral Monte Carlo quantum annealing algorithm; (2) Initialize P quadratic unconstrained binary optimization matrices QUBO to form the original matrix set E and store them in P BRAM1s in the FPGA. (3) Coupling coefficient J g , read the original matrix set E from BRAM1, generate random numbers through FPGA, use the random numbers to perform random operations on the original matrix set E to update P new QUBO matrices, form a new matrix set E′, and store it in another P BRAM2 in FPGA; (4) Calculate the potential energy difference ΔH between the original matrix set E and the new matrix set E′ respectively p , kinetic energy difference ΔH k , the difference of the classical Hamiltonian ΔH c ; (5) For the potential energy difference ΔH p and the classical Hamiltonian difference ΔH c Compare, if one of the differences is less than 0, or the classical Hamiltonian difference ΔH c If the Monte Carlo criterion is met, the corresponding matrix in the original matrix set E is replaced by the corresponding matrix in the new matrix set E′, and the matrix with the smallest potential energy term in the updated original matrix set E is selected to obtain the optimal solution s of the current Monte Carlo cycle best ; (6) Update the coupling coefficient J g , the number of iterations M minus 1, to determine whether the current number of iterations M is greater than 0: If yes, return to step (3); Otherwise, the FPGA stops working and obtains the optimal path result, which is the matrix s of the current iteration. best .

2. The method according to claim 1, characterized in that In step (2), P quadratic unconstrained binary optimization matrices QUBO are initialized to form the original matrix set E, including the following: (2a) The original matrix has N rows and N columns. All elements in the matrix have only two values: 0 and 1. The rows represent the node numbers in the path planning, and the columns represent the time points. (2b) In the current row, randomly select an integer a between 0 and N and assign the element in the ath column of the current row to 1; (2c) Repeat (2b) until all rows are selected to obtain a separate original matrix s; (2d) Repeat (2a) to (2c) until all P matrices are initialized to form the original matrix set E; E={s 1 ,s 2 ,…s m …,s P }, where s m The mth original matrix s in the original matrix set E, where m is an integer between 1 and P; each matrix is ​​a randomly arranged matrix of 0s and 1s.

3. The method according to claim 1, characterized in that In step (3), the random number generated by FPGA is generated by using a true random seed and an m sequence to generate a random integer between 0 and N that satisfies a uniform distribution. The implementation includes the following: (3a) Design an oscillator consisting of three series-connected NOT gates on the FPGA backplane; (3b) interconnecting two oscillators in pairs to form an oscillation ring. At the sampling clock moment, the output of the first oscillator a is connected to the input of the second oscillator b, and the output of the second oscillator b is simultaneously connected to the input of the first oscillator a. (3c) Arrange 19 oscillator rings in parallel in the FPGA, perform XOR operation on all output results, and output single-bit data; (3d) Use 32 registers in the FPGA to collect single-bit data for 32 clock cycles to obtain an initial random number seed as the initial state of the m-sequence; (3e) Predict the value of the m sequence after log2(N) clock cycles; (3f) Truncate the lower log2(N) bits of the newly generated m sequence to obtain a random number that meets the requirements.

4. The method according to claim 1, characterized in that Step (3) performs random operations on the original matrix set S using random numbers, including: (3g) Take two different random row numbers c and f in sequence, and determine whether the second random number d is an odd number: If so, all columns between c and f are flipped symmetrically; Otherwise, the selected rows c and f are swapped.

5. The method according to claim 1, characterized in that The new matrix set E′ in step (3) is E′={s ′1 ,s ′2 ,…s ′m …,s ′P }, where s ′m Represents the mth new matrix s′ in the new matrix set E′, where m is an integer from 1 to P; each matrix is ​​a randomly arranged matrix of 0s and 1s with N rows and N columns.

6. The method according to claim 1, characterized in that In step (4), the difference ΔH between the potential energy term of the original matrix set S and the new matrix set S′ is calculated. p ,include: (4a) Calculate the potential energy term of the original matrix Among them, d i,j It is the element with horizontal coordinate i and vertical coordinate j in the position matrix d, which represents the distance between node i and node j; Indicates m The element with horizontal coordinate i and vertical coordinate k; s mT Yes m The transition matrix obtained by shifting left by 1 bit is Denotes the transition matrix s mT The element with horizontal coordinate k and vertical coordinate j; Indicates m The element with horizontal coordinate i and vertical coordinate j; W max Represents the element with the largest value in the position matrix d; (4b) Calculate the potential energy term of the new matrix Among them, s ′m is the mth new matrix in the new matrix set E′ obtained after random operation; Indicates ′m The element with horizontal coordinate i and vertical coordinate k; s′ mT s′ m The new transition matrix obtained by shifting left by 1 bit; for The element with horizontal coordinate k and vertical coordinate j; Indicates ′m The element with horizontal coordinate i and vertical coordinate j; W max Represents the element with the largest value in the position matrix d. (4c) According to and The potential energy difference ΔH is obtained p :

7. The method according to claim 1, characterized in that In step (4), the difference ΔH between the kinetic energy terms of the original matrix set S and the new matrix set S′ is calculated: k ,include: (4d) Calculate the kinetic energy term of the original matrix in, Represents the element with horizontal coordinate i and vertical coordinate j in the m-1th matrix in the original matrix set E; Represents the element with horizontal coordinate i and vertical coordinate j in the mth matrix in the original matrix set E; Represents the element with horizontal coordinate i and vertical coordinate j in the m+1th matrix in the original matrix set E; (4e) Calculate the kinetic energy term of the new matrix in, represents the element with horizontal coordinate i and vertical coordinate j in the m-1th matrix in the new matrix set E′; represents the element with horizontal coordinate i and vertical coordinate j in the mth matrix in the new matrix set E′; represents the element with horizontal coordinate i and vertical coordinate j in the m+1th matrix in the new matrix set E′; (4f) According to and The kinetic energy difference ΔH is obtained k :

8. The method according to claim 1, characterized in that: The coupling coefficient J is calculated in step (3) g : I g =sin(-T / 2)*log(tanh(Γ / PT)); Calculate the classical Hamiltonian difference ΔH in step (4) c , the formula is as follows: ΔH c =ΔH p / P-J g ×ΔH k 。 9. A Monte Carlo path integral-simulated quantum annealing system based on FPGA, comprising: Bram storage space, used to store the original matrix set and the new matrix set; Bram space access module, used to control the read and write addresses of Bram storage space; A random number transmitting module is used to generate a random number between 0 and N; Random operation module, which performs random operations on the original matrix according to random numbers to obtain a new matrix; Potential energy term calculation module, used to solve the potential energy terms of the original matrix and the new matrix; Kinetic energy term calculation module, used to solve the kinetic energy terms of the original matrix and the new matrix; The judgment update module is used to calculate the interpolation of potential energy terms, the difference of kinetic energy terms, and the difference of classical Hamiltonian, and judge whether to replace the old matrix with the new matrix based on these results. The iteration counting module is used to reduce the number of iterations M by 1 after the update module is executed, and to determine whether the number of iterations M is greater than 0: if so, return to the random operation module; Otherwise: FPGA stops working.

10. An electronic device based on FPGA, comprising a memory, a reconfigurable logic operation unit and a computer program that can be run on the FPGA, characterized in that: When the FPGA executes the program, the Monte Carlo path integral-simulated quantum annealing method according to any one of claims 1 to 9 is implemented.

Images