Tokamak plasma current density reconstruction system based on a coherent ising machine
By combining electromagnetic diagnostic data acquisition, Bayesian posterior modeling, and coherent Ising machine quantum solution, the large-scale tokamak plasma current density reconstruction problem is decomposed into multiple sub-problems, solving the reconstruction problem under the hardware limitations of the coherent Ising machine and achieving efficient and accurate plasma current density reconstruction.
Patent Information
- Authority / Receiving Office
- CN · China
- Patent Type
- Applications(China)
- Current Assignee / Owner
- INST OF ENERGY HEFEI COMPREHENSIVE NAT SCI CENT (ANHUI ENERGY LAB)
- Filing Date
- 2026-06-04
- Publication Date
- 2026-07-03
AI Technical Summary
Existing technologies struggle to efficiently solve the large-scale continuous quadratic optimization problem of tokamak plasma current density reconstruction under the limited bit capacity of the coherent Ising machine hardware. Furthermore, existing methods have high computational complexity and cannot meet the requirements of real-time control systems.
A tokamak plasma current density reconstruction system based on the coherent Ising machine is adopted. By combining electromagnetic diagnostic data acquisition, Bayesian posterior modeling, spatial block processing, block coordinate descent iterative control, and the coherent Ising machine quantum solver, the large-scale optimization problem is decomposed into multiple sub-problems and solved in parallel to achieve efficient reconstruction.
It achieves high-precision and rapid reconstruction of plasma current density profiles under the constraint of existing quantum hardware bit capacity, and quantifies reconstruction uncertainty, which significantly improves reconstruction speed and accuracy and meets the real-time control requirements of tokamak devices.
Smart Images
Figure CN122334540A_ABST
Abstract
Description
Technical Field
[0001] This invention relates to the interdisciplinary field of quantum computing and nuclear fusion plasma profile reconstruction, and particularly to a tokamak plasma current density reconstruction system based on the coherent Ising machine. Background Technology
[0002] The EAST (Expert Superconducting Tokamak) is my country's independently developed world's first fully superconducting tokamak device. Its core objective is to achieve high-parameter steady-state plasma operation, laying the engineering foundation for future fusion reactors. In a tokamak device, the plasma current density profile is one of the most important control parameters, directly determining the plasma configuration, stability, and confinement performance. Accurately reconstructing the plasma current density profile is fundamental to achieving plasma configuration control and studying magnetohydrodynamic instabilities. However, this profile cannot be directly measured and must be reconstructed through inverse inference from external electromagnetic measurement data (magnetic probes, flux loops, Rochelle coils, etc.). This is an underdetermined, ill-conditioned, large-dimensional inverse problem, representing a core technical challenge in the field of nuclear fusion plasma control.
[0003] Currently, the mainstream methods for reconstructing plasma current density profiles mainly include the following two:
[0004] (1) EFIT Method: EFIT (Equilibrium Fitting) is the most widely used classical method in the field of tokamak equilibrium reconstruction, proposed by Lao et al. in 1985. This method is based on the magnetohydrodynamic equilibrium equation (Grad-Shafranov equation). After discretization on the poloidal section, it uses the least squares fitting method to invert the plasma pressure profile and poloidal current function using external magnetic measurement data, and then reconstructs the magnetic flux surface shape and current density distribution. EFIT has been widely adopted by major tokamak devices in the world such as DIII-D, JET, KSTAR, and EAST, and there is a real-time version (RT-EFIT) for closed-loop configuration control. However, the EFIT method has the following shortcomings: it relies on parametric profile assumptions (such as polynomial or spline functions), lacks flexibility, and is difficult to accurately capture non-smooth or non-monotonic profile features that may occur in the experiment; it is sensitive to initial values and easily gets trapped in local extrema; and it cannot quantitatively give the reconstruction uncertainty.
[0005] (2) Classical numerical method based on Bayesian inference: This method models the plasma current density reconstruction problem as a Bayesian posterior probability maximization problem. By constructing prior probabilities (describing the spatial smoothness and physical rationality of the plasma current density) and likelihood probabilities (fitting external electromagnetic measurement data), it utilizes maximum a posteriori estimation to simultaneously consider data fitting accuracy and physical rationality, and can provide a quantification of statistical uncertainty. The core of this method is a large-scale quadratic optimization problem: ,in Let be the plasma current density vector, Q be the posterior precision matrix (positive definite), and c be the linear coefficient vector. The computational complexity of the classical numerical solution is O(n log n). (Cholesky decomposition) is far from meeting the reconstruction speed requirements of the real-time control system of the tokamak device.
[0006] In recent years, quantum computing technology has made rapid progress. The Coherent Ising Machine (CIM) is a dedicated quantum simulation hardware based on an optical parametric oscillator, capable of efficiently solving Ising models or equivalent quadratic unconstrained binary optimization (QUBO) problems. The main advantages of CIM are its strong parallelism and outstanding ability to solve densely coupled problems. However, directly transforming this large-scale continuous quadratic optimization problem into a QUBO problem (requiring multiple bits per variable after binarization, resulting in a total bit count exceeding 3000) would far exceed the bit capacity limit of existing CIM hardware, making quantum acceleration impractical in engineering. Therefore, this application proposes a tokamak plasma current density reconstruction system based on the Coherent Ising Machine. Summary of the Invention
[0007] The purpose of this invention is to address the large-scale continuous quadratic optimization problem in the prior art of reconstructing plasma current density under the constraint of limited hardware bit capacity of the coherent Ising machine, and to propose a tokamak plasma current density reconstruction system based on the coherent Ising machine.
[0008] The technical solution of this invention: a tokamak plasma current density reconstruction system based on a coherent Ising machine, comprising:
[0009] The electromagnetic diagnostic data acquisition unit is used to acquire external electromagnetic measurement data of the tokamak device in real time and form a measurement vector.
[0010] The Bayesian posterior modeling processor, connected to the electromagnetic diagnostic data acquisition unit, is used to construct the posterior accuracy matrix and linear term coefficient vector based on the conditional autoregressive prior model and the electromagnetic measurement likelihood probability, transforming the plasma current density reconstruction problem into a standard quadratic optimization problem.
[0011] The spatial block processor is used to acquire the coordinate information of active grid points in the plasma poloidal section, divide the active grid points into multiple sub-blocks according to a preset spatial block strategy, and output the variable index set and corresponding sub-matrix of each sub-block.
[0012] The block coordinate descent iterative controller is connected to the spatial block processor and the Bayesian posterior modeling processor respectively. It is used to execute the block coordinate descent iterative process with momentum acceleration. In each iteration, a corresponding quadratic unconstrained binary optimization subproblem is generated for each sub-block.
[0013] A coherent Ising machine quantum solver, connected to a block coordinate descent iterative controller, is used to solve in parallel the various quadratic unconstrained binary optimization subproblems generated by the block coordinate descent iterative controller and output the spin configuration with the lowest energy.
[0014] The results are output to the visualization unit, which is connected to the block coordinate descent iterative controller to receive the final solved current density vector and map it to generate a current density profile and reconstruct the uncertainty estimate.
[0015] Optionally, the Bayesian posterior modeling processor is specifically used for:
[0016] A conditional autoregressive prior model is adopted, the adjacency matrix is defined based on the adjacency relationship of active grid points, and a prior accuracy matrix is constructed.
[0017] A linear positive model of the response matrix of the electromagnetic diagnostic system is constructed based on the Biot-Savart law, and a likelihood probability model is constructed assuming that the measurement noise follows a zero-mean Gaussian distribution.
[0018] By combining the conditional autoregressive prior with Gaussian likelihood, the posterior probability density is derived to be a Gaussian distribution using Bayes' theorem.
[0019] The problem of maximizing the posterior probability is equivalently transformed into minimizing the standard quadratic objective function, resulting in the posterior accuracy matrix and the linear term coefficient vector.
[0020] Optionally, the spatial partitioning processor is specifically used for:
[0021] Extract the active grid points from the discrete grid of the poloidal section of the tokamak device to solve for the current density;
[0022] A uniform block partitioning strategy is adopted, which divides the coordinate space into several regions according to the preset block width and block height, and assigns active grid points located in the same region to the same spatial sub-block;
[0023] For sub-blocks with fewer than a preset minimum block threshold, they are merged with adjacent normal-sized sub-blocks to ensure that the size of each sub-block meets the solution efficiency requirements of the coherent Ising machine quantum solver.
[0024] Optionally, the partitioning strategy employed by the spatial partitioning processor can be replaced by any of the following strategies:
[0025] The striping block strategy divides active grid points into blocks of equal length in the order of the global index.
[0026] The graph partitioning strategy is to divide the graph composed of active grid points into several subgraphs by using the absolute values of the elements of the posterior precision matrix as edge weights, so as to minimize the total weight of the cut edges.
[0027] An adaptive block partitioning strategy is adopted, which dynamically adjusts the block size based on the number of available bits in the coherent Ising machine quantum solver.
[0028] Optionally, the block coordinate descent iterative controller includes an initialization module, a main update module, a sub-iteration refinement module, a momentum extrapolation module, and a convergence judgment module, wherein:
[0029] The main update module is used to extract the local objective function matrix and coupling vector of each sub-block, and encode them into a quadratic unconstrained binary optimization sub-problem and submit it to the coherent Ising machine quantum solver for solving.
[0030] The sub-iterative refinement module is used to fix the integer part of the solution and use the fractional part as a variable to reconstruct and solve the quadratic unconstrained binary optimization sub-problem with higher bit precision based on the current master update solution, so as to improve the accuracy of the solution.
[0031] The momentum extrapolation module is used to calculate the momentum coefficient based on the momentum strategy of the fast iterative shrinkage threshold algorithm, and to apply extrapolation to the current solution to accelerate convergence.
[0032] Optionally, the block coordinate descent iterative controller further includes a quadratic unconstrained binary optimization encoding module, used to perform the following operations before submitting the subproblem to the coherent Ising machine quantum solver:
[0033] Estimate the range of values for each continuous variable based on the eigenvalues of the posterior precision matrix;
[0034] Each continuous variable is encoded using an integer quantum variable, and represented as a weighted sum of multiple binary bits;
[0035] Substituting the binary code into the local objective function of the sub-block, we can expand it to obtain a quadratic unconstrained binary optimization form.
[0036] The quadratic unconstrained binary optimization problem is transformed into a standard Ising model through variable transformation, and then solved by the coherent Ising machine quantum solver.
[0037] Optionally, the conditional autoregressive prior model used by the Bayesian posterior modeling processor can be replaced by any of the following models:
[0038] The advanced squared exponential prior model, which uses a Gaussian process kernel function to describe the spatial correlation of current density;
[0039] Laplace prior model, used to induce sparsity of solutions;
[0040] The mixture of Gaussian process prior model, which is composed of multiple Gaussian process priors, is used to capture multi-scale spatial correlations.
[0041] Optionally, the coherent Ising machine quantum solver is a dedicated quantum simulation hardware based on an optical parametric oscillator; when the coherent Ising machine quantum solver is unavailable, the block coordinate descent iterative controller seamlessly switches to a classical solver, which is one or more combinations of fast simulated annealing algorithm, D-Wave quantum annealing machine, branch and bound method or tabu search algorithm.
[0042] Optionally, the tokamak device is an EAST device, and the external electromagnetic measurement data acquired by the electromagnetic diagnostic data acquisition unit includes signals from a magnetic probe, a flux loop, and a Rochelle coil; the system decomposes a large-scale quadratic optimization problem with an overall dimension of hundreds of dimensions into multiple sub-problems with dimensions of tens of dimensions, so that the number of quadratic unconstrained binary optimization bits after encoding each sub-problem does not exceed the upper limit of the hardware bit capacity of the coherent Ising machine quantum solver.
[0043] Compared with the prior art, this application includes at least one of the following beneficial technical effects:
[0044] This invention decomposes a large-scale quadratic optimization problem of approximately 420 dimensions into multiple subproblems of approximately 12 to 20 dimensions through the collaborative work of a spatial block processor and a block coordinate descent iterative controller. Each subproblem, after being encoded by a second unconstrained binary optimization, does not exceed the capacity limit of approximately 550 bits of the existing coherent Ising machine hardware. This enables quantum computing resources to be effectively adapted to the large-scale physical inverse problem of plasma current density reconstruction for the first time.
[0045] This invention constructs a posterior accuracy matrix through a Bayesian posterior modeling processor and combines it with the FISTA momentum acceleration mechanism to achieve a convergence rate of O(1 / k²) for the block coordinate descent algorithm. Compared with the convergence rate of O(1 / k) for ordinary block coordinate descent, the number of iterations can be reduced by more than 50%, and the maximum relative error of plasma current density reconstruction can be controlled within 10%.
[0046] After the main iteration of block coordinate descent is completed, this invention fixes the integer part of the current solution, re-encodes the fractional part with higher bit precision, and submits it to the coherent Ising machine for solution. This significantly improves the solution quality in the later stages of convergence without increasing the complexity of the main iteration.
[0047] When the coherent Ising machine quantum solver is unavailable, the block coordinate descent iterative controller can seamlessly switch to classical solvers such as fast simulated annealing, ensuring stable operation of the system in real experimental environments.
[0048] In summary, this invention organically integrates Bayesian inference modeling, spatial block decomposition, block coordinate descent iterative control, and coherent Ising machine quantum solution, achieving high-precision and rapid reconstruction of tokamak plasma current density profiles with uncertainty quantification under the constraint of existing quantum hardware bit capacity. Attached Figure Description
[0049] Figure 1 This is a system architecture diagram for reconstructing the plasma current density of a tokamak based on a coherent Ising machine;
[0050] Figure 2 This is a schematic diagram of a 33×33 active point spatial block. Detailed Implementation
[0051] The following specific examples illustrate the implementation of the present invention. Those skilled in the art can easily understand other advantages and effects of the present invention from the content disclosed in this specification. The present invention can also be implemented or applied through other different specific embodiments, and various details in this specification can also be modified or changed based on different viewpoints and applications without departing from the spirit of the present invention. It should be noted that, unless otherwise specified, the following embodiments and features described therein can be combined with each other.
[0052] Example 1: The tokamak plasma current density reconstruction system based on the coherent Ising machine proposed in this invention will be described in detail below.
[0053] I. System Overall Architecture
[0054] like Figure 1 This invention provides a tokamak plasma current density reconstruction system based on a coherent Ising machine. The system consists of six functional modules, which are connected sequentially via a data bus:
[0055] (1) Electromagnetic diagnostic data acquisition unit: responsible for real-time acquisition of external electromagnetic measurement data of the EAST tokamak device, including the output signals of magnetic probes (Pickup Coils), flux loops, and Rogowski coils, forming... 3D measurement vector And complete analog-to-digital conversion and noise filtering preprocessing. ,in Represents current Electromagnetic measurement The response coefficients between Indicates other current sources to the measurement vector Contribution To represent measurement error, assume It follows a Gaussian distribution.
[0056] (2) Bayesian posterior modeling processor: Receives the output of the electromagnetic diagnostic data acquisition unit, and transforms the plasma current density reconstruction problem into a standard quadratic optimization problem based on the CAR prior model and the electromagnetic measurement likelihood probability. .
[0057] in Let be the posterior precision matrix, and be a positive definite matrix. Represents the vector of linear term coefficients. This is a constant term.
[0058] (3) Spatial Block Processor: Receives the coordinate information of active grid points, groups the active points into several sub-blocks according to the spatial block strategy, and outputs the variable index set and corresponding sub-matrix of each sub-block. , To the block coordinate descent iterative controller, where yes Matrix in sub-block The inner diagonal blocks describe the coupling relationships between variables within the block. yes Sub-blocks in the matrix The off-diagonal blocks between the other variable blocks will condense cross-block coupling into effective linear terms.
[0059] (4) Coherent Ising Machine (CIM) Quantum Solver: The core hardware unit of this system, based on an optical parametric oscillator, used to solve the QUBO / Ising subproblems in parallel. It receives the Ising matrix encoded by the block coordinate descent iterative controller and outputs the spin configuration with minimum energy.
[0060] (5) Block coordinate descent iterative controller: coordinates the data flow and control logic between modules, executes the BCD main iterative process with FISTA momentum acceleration, manages sub-iteration refinement and convergence judgment, and drives the CIM quantum solver to solve the sub-block QUBO problem in a loop when convergence is not achieved.
[0061] (6) Result Output and Visualization Unit: Receives the final current density vector output by the block coordinate descent iterative controller. This is mapped back to the profile grid to form a current density profile, and the reconstruction uncertainty estimate is output.
[0062] II. Working Principle of Bayesian Posterior Modeling Processor
[0063] The Bayesian posterior modeling processor internally implements the following computational logic to transform electromagnetic measurement data into a standard quadratic optimization objective:
[0064] (1) Prior probability Modeling: This system employs a Conditional Autoregressive (CAR) prior model to describe the spatial smoothness of the plasma current density. The CAR prior is based on grid adjacency relationships, defining an adjacency matrix. (If grid points) and Adjacent (otherwise it is 0), precision matrix The following formula is defined:
[0065]
[0066]
[0067] in, Variance is the distance between the quantized value and the mean. For spatially relevant parameters (take) To ensure positive definiteness). Let be the identity matrix. Based on the CAR prior, the prior probability of the current vector. The density follows a Gaussian distribution.
[0068]
[0069] in, For the current vector dimension, for The size of the prior covariance matrix.
[0070] (2) Likelihood probability modeling: Measurement data vector of EAST electromagnetic diagnostic system It can be represented as:
[0071]
[0072] in, This is the output signal of the magnetic probe. This is the output signal of the flux loop. The output signal of the Roco coil;
[0073] According to the Biot-Savart law, the magnetic field generated by a current element is proportional to the current. All diagnostic response coefficients are integrated into the response matrix. In the process, a positive model is obtained. :
[0074]
[0075] in, The response matrix describes the linear relationship between the current beam and the electromagnetic measurement. This is a vector of constant terms (such as the contribution of the PF coil current to the diagnostic signal). For measuring noise, it follows a zero-mean Gaussian distribution. Therefore, the likelihood probability density is:
[0076]
[0077] in, To measure the covariance matrix, it is usually a diagonal matrix, with elements representing the variance of each measurement channel. Representation matrix Dimensions This represents the transpose of a matrix.
[0078] (3) Posterior probability and objective function: the prior probability Likelihood probability Substituting into Bayes' theorem, we obtain the posterior probability density:
[0079]
[0080] in, Let be the posterior probability, i.e., the objective function.
[0081] Maximize posterior probability This is equivalent to minimizing the negative of the logarithmic posterior function, i.e., minimizing the following quadratic objective function:
[0082]
[0083] in Let be the posterior precision matrix, and be a positive definite matrix. Represents the vector of linear term coefficients. This is a constant term. Because... and All are positive definite matrices. It is also a positive definite matrix.
[0084] III. Working Principle of Spatial Block Processor
[0085] The spatial block processor decomposes approximately 420 active grid points into several sub-blocks that can be processed independently by the CIM hardware, as follows:
[0086] (1) Active point extraction: The poloidal section of the EAST tokamak was discretized into a uniform 33×33 grid. Approximately 669 grid points located outside the plasma boundary or confinement area were considered frozen points with a fixed current density of zero. The remaining approximately 420 grid points were considered active points, with their current density being the variable to be solved. The active grid points were extracted by reading their coordinates. And the corresponding index, construct the active point coordinate sequence.
[0087] (2) Block division: according to block width and block height The coordinate space is uniformly divided, and each active point is assigned to:
[0088] same Active points are grouped into the same spatial sub-block. Operators The floor operator is mathematically defined as: for any real number... , Indicates less than or equal to The largest integer. Typical parameters are set to... , At this point, the active region is divided into several sub-blocks of similar size, each containing approximately 12 to 20 active points.
[0089] (3) Merging excessively small blocks: For sub-blocks with fewer active points than the minimum block size threshold min_block_size (default value 5), prioritize merging them with normally sized 4-adjacent (top, bottom, left, right) sub-blocks to ensure that each sub-block is of reasonable size and avoid excessively small subproblems that lead to wasted CIM resources or solution degradation. After merging, according to... The spatial order renumbers all sub-blocks, such as Figure 2 The diagram shows a 33×33 active point spatial block.
[0090] IV. Working Process of the Block Coordinate Descent Iterative Controller
[0091] The block coordinate descent iterative controller is the scheduling core of this system, executing the BCD iterative process with FISTA momentum acceleration. The complete workflow is as follows:
[0092] Phase 1, Initialization: Initial Solution It can be generated by a zero vector, a random vector, or a fast simulated annealing (FSA) hot start; initialize momentum variables. Time parameters .
[0093] Phase 2, Sub-block Master Update: For each sub-block Execute in sequence: (a) Spatial block processor extracts sub-blocks Local objective function matrix (a) The coupling vector; (b) The QUBO encoding module converts the continuous variables into binary bit representations and constructs the Ising matrix; (c) The Ising matrix is submitted to the CIM quantum solver for parallel solution, and the spin combination with the lowest energy is selected; (d) The continuous variable values are decoded and recovered. If the objective function value decreases, then an update is accepted.
[0094] Phase 3, Sub-iteration Refinement: For each sub-block, after Phase 2 is completed, ... Represents the integer part of the current density solution. ,by This indicates that the decimal part is treated as a new variable. , ,use Replace the variables in the original QUBO problem Reconstruct the QUBO problem and resubmit it for CIM solution. If the objective function value decreases after refinement, then accept the refined solution.
[0095] Phase 4, FISTA Momentum Extrapolation: Calculating Momentum Parameters momentum coefficient Apply extrapolation to the current solution If the objective function value increases after extrapolation, then reset. , .in It is the first The momentum-time parameter of the next iteration controls the extrapolation step size; It is the momentum coefficient, used to calculate extrapolated quantities; It is an auxiliary variable derived from momentum extrapolation; It is the first The current density solution obtained in the second iteration.
[0096] Phase 5, Convergence Criterion: If the objective function decreases by a certain amount... Or reach the maximum number of iterations, where It is the first The objective function value of the next iteration. It is the convergence threshold of the objective function, preset to 1. The iterative controller outputs the current optimal solution to the result output unit; otherwise, it triggers a de-stressing action and returns to stage two to continue iterating.
[0097] V. CIM Quantum Solver and QUBO Encoding Module
[0098] After receiving the Ising matrix output by the QUBO encoding module, the CIM quantum solver searches in parallel for the lowest-energy spin configuration using an optical parametric oscillation network. The working steps of the QUBO encoding module are as follows:
[0099] Step 1, Variable range estimation: based on the posterior accuracy matrix The smallest eigenvalue with posterior mean Estimate the upper bound of the norm of the optimal solution. Based on this, the range of values for each variable is determined. .in For the posterior accuracy matrix The smallest eigenvalue; This is the posterior mean vector; The optimal current density solution; A constant vector related to the observed data; and The first The lower and upper bounds of each variable.
[0100] Step 2, Binary Encoding of Continuous Variables: For each continuous variable within a sub-block... Encode integer quantum variables: ,in Indicates the first A continuous current density variable, The number of bits (typically 8). Scaling factor It is a binary variable.
[0101] Step 3, QUBO matrix construction: Substitute the binary code into the local objective function Expand to obtain the QUBO form ,in, for Vie Coupling matrix. For sub-blocks The vector of current density variables; For sub-blocks The precision matrix (i.e.) ); For sub-blocks The coupling matrix between it and the remaining sub-blocks; To remove child blocks The current density variable vector of the remaining sub-blocks; For sub-blocks The linear term coefficient vector; For sub-blocks The number of variables in the text; For sub-blocks A binary encoded variable vector, This represents the transpose of a matrix.
[0102] Step 4, QUBO to Ising: through variable transformation ( ), converted into standard Ising form If the size of the Ising matrix exceeds the CIM hardware bit capacity limit. Then, adaptive rounding is performed on the matrix precision, where It is the energy function of the Ising model, representing the spin configuration. The corresponding system energy; For the first Each Ising spin variable takes values from _____. , indicating the first The spin state of each qubit is represented by -1 (spin down) and +1 (spin up). For the first One Ising spin variable; It is a QUBO binary variable, taking values from 1 to 2. ; In the Ising model, the first Local field coefficients (linear terms) of each spin; In the Ising model, the first With the The coupling coefficient between the spins; This represents the upper limit of the bit capacity of the CIM hardware.
[0103] Step 5, CIM Solving and Decoding: The CIM quantum solver runs in parallel on the same Ising matrix. Second-rate, This represents the number of times the CIM quantum solver performs parallel operations on the same Ising matrix, i.e., the number of samplings, used to statistically obtain the spin configuration with the lowest energy. The spin configuration with the lowest energy is then selected, and an inverse transformation is performed to recover the continuous variable values.
[0104] VI. Alternative Solutions
[0105] 1. Alternative solutions to the block-based strategy
[0106] The spatial block processor of this system discloses a block partitioning strategy based on spatial coordinates, but the following alternative partitioning strategies can also be used:
[0107] (1) Strip method (sequential block division): Divide according to the global sequence number, and each segment is continuous. Each active point forms a block, requiring no coordinate information, and is applicable to any grid topology.
[0108] (2) Graph partitioning method: based on the a posteriori precision matrix absolute value of elements To determine edge weights, the variable node graph is divided into several subgraphs to minimize the total weight of the cut edges. This approach is suitable for... The case where the matrix has good sparsity.
[0109] (3) Adaptive block division: The block size is dynamically adjusted according to the number of available bits in the CIM hardware. When there are enough available bits, the block size is increased to reduce the number of iterations.
[0110] 2. Alternative solutions to prior models
[0111] The Bayesian posterior modeling processor in this system uses a CAR (Conditional Autoregressive) prior model, but the following alternative prior models can also be used:
[0112] (1) ASE Prior (Advanced Squared Exponential Prior): The spatial correlation of current density is described using a Gaussian process kernel function. The ASE prior covariance matrix is defined as follows: ,in The first of the prior covariance matrices of ASE The element represents the element. The grid point and the first The covariance between the current density values at each grid point; This is a signal variance hyperparameter that controls the overall amplitude of current density fluctuations. For the first The grid point and the first Euclidean distance (spatial distance) between grid points; The correlation length hyperparameter controls the decay rate of the spatial correlation of current density. A larger value indicates a slower decay of spatial correlation. This prior has a better condition number, making it suitable for handling scenarios where the current density distribution has multi-scale characteristics.
[0113] (2) Laplace prior: Using the Laplace distribution as a prior can induce the sparsity of the solution, which is suitable for scenarios where there are local peaks in plasma current density.
[0114] (3) Mixed Gaussian process priors: By superimposing multiple Gaussian process priors, multi-scale spatial correlations can be captured simultaneously, improving the reconstruction accuracy of complex current profiles.
[0115] 3. Alternatives to the CIM quantum solver
[0116] This system uses CIM as the primary quantum solver, and Fast Simulated Annealing (FSA) is used as a classical alternative when CIM is unavailable. In addition, the following solvers can also be embedded in this system:
[0117] (1) D-Wave quantum annealing machine: The Ising matrix is embedded in the Pegasus graph topology of D-Wave and then solved. It is suitable for subproblems of sparse coupled matrices.
[0118] (2) Classical Branch and Bound Method: For small-scale subproblems , Indicates a sub-block The number of variables in the sub-block refers to the total number of current density variables (or their binary encoded forms) that need to be optimized within that sub-block. The branch and bound method can be used to solve this problem precisely.
[0119] (3) Tabu Search: A metaheuristic classical optimization algorithm suitable for solving large-scale QUBO problems. It can be used as an alternative to FSA when CIM is unavailable.
[0120] 4. Sub-iteration refinement alternatives
[0121] This system employs a sub-iteration refinement strategy (fixing the integer part after the main iteration and performing high-precision QUBO encoding to solve the fractional part) to improve the accuracy of the solution. Alternative solutions include:
[0122] (1) Direct high-precision solution: directly using higher bit precision (e.g. (or higher) QUBO encoding is used to solve all variables without distinguishing between main iterations and sub-items. This method is simple and direct but requires a large number of bits, and is only suitable for scenarios where CIM hardware has sufficient bit capacity.
[0123] (2) Step-by-step refinement: First, refine with low precision (coarse bit count, such as...) The method obtains a coarse-grained solution, and then iteratively refines the precision by increasing the precision by 1 bit each time, using this solution as the center, until the target precision is reached. This method is more stable than sub-iterative refinement but requires more iterations, making it suitable for scenarios with ample hardware bit resources.
[0124] It is worth noting that this invention, through the collaborative work of a spatial block processor and a block coordinate descent iterative controller, decomposes a 420-dimensional large-scale quadratic optimization problem into several QUBO subproblems of approximately 12 to 20 dimensions. Each subproblem, after QUBO encoding, has a size not exceeding 550 bits, fully satisfying the constraints of existing CIM hardware. Compared to the direct global QUBO scheme, which requires more than 3000 bits, this system can be implemented on existing CIM hardware, enabling quantum computing resources to be effectively adapted to the large-scale physical inverse problem of plasma current density reconstruction for the first time.
[0125] The Bayesian posterior modeling processor, based on CAR prior and quantum BCD iterative optimization with FISTA momentum acceleration, can reduce the maximum relative error of plasma current density reconstruction to within 10%. Compared to a regular BCD The number of iterations is reduced by more than 50%. The sub-iteration refinement mechanism further reduces the objective function value, improving solution quality without increasing the complexity of the main iteration.
[0126] This system utilizes a sub-iteration refinement module within the block coordinate descent iterative controller. Building upon the fixed integer portion of the main iteration, it performs higher-precision QUBO encoding on the fractional part before submitting it to the CIM solver. This mechanism significantly improves solution accuracy in the later stages of convergence without increasing the complexity of the main iteration. Furthermore, the CIM / FSA dual solver support (seamlessly switching to fast simulated annealing when CIM hardware is unavailable) ensures reliable operation of the system in the EAST real-world experimental environment and can be extended to other tokamak devices such as HL-2A and CFETR.
[0127] The above specific embodiments are merely several optional embodiments of the present invention. Based on the technical solutions of the present invention and the relevant teachings of the above embodiments, those skilled in the art can make various alternative improvements and combinations to the above specific embodiments.
Claims
1. A tokamak plasma current density reconstruction system based on a coherent Ising machine, characterized in that, include: The electromagnetic diagnostic data acquisition unit is used to acquire external electromagnetic measurement data of the tokamak device in real time and form a measurement vector. The Bayesian posterior modeling processor, connected to the electromagnetic diagnostic data acquisition unit, is used to construct the posterior accuracy matrix and linear term coefficient vector based on the conditional autoregressive prior model and the electromagnetic measurement likelihood probability, transforming the plasma current density reconstruction problem into a standard quadratic optimization problem. The spatial block processor is used to acquire the coordinate information of active grid points in the plasma poloidal section, divide the active grid points into multiple sub-blocks according to a preset spatial block strategy, and output the variable index set and corresponding sub-matrix of each sub-block. The block coordinate descent iterative controller is connected to the spatial block processor and the Bayesian posterior modeling processor respectively. It is used to execute the block coordinate descent iterative process with momentum acceleration. In each iteration, a corresponding quadratic unconstrained binary optimization subproblem is generated for each sub-block. A coherent Ising machine quantum solver, connected to a block coordinate descent iterative controller, is used to solve in parallel the various quadratic unconstrained binary optimization subproblems generated by the block coordinate descent iterative controller and output the spin configuration with the lowest energy. The results are output to the visualization unit, which is connected to the block coordinate descent iterative controller to receive the final solved current density vector and map it to generate a current density profile and reconstruct the uncertainty estimate.
2. The tokamak plasma current density reconstruction system based on the coherent Ising machine according to claim 1, characterized in that, The Bayesian posterior modeling processor is specifically used for: A conditional autoregressive prior model is adopted, the adjacency matrix is defined based on the adjacency relationship of active grid points, and a prior accuracy matrix is constructed. A linear positive model of the response matrix of the electromagnetic diagnostic system is constructed based on the Biot-Savart law, and a likelihood probability model is constructed assuming that the measurement noise follows a zero-mean Gaussian distribution. By combining the conditional autoregressive prior with Gaussian likelihood, the posterior probability density is derived to be a Gaussian distribution using Bayes' theorem. The problem of maximizing the posterior probability is equivalently transformed into minimizing the standard quadratic objective function, resulting in the posterior accuracy matrix and the linear term coefficient vector.
3. The tokamak plasma current density reconstruction system based on the coherent Ising machine according to claim 1, characterized in that, The spatial block processor is specifically used for: Extract the active grid points from the discrete grid of the poloidal section of the tokamak device to solve for the current density; A uniform block partitioning strategy is adopted, which divides the coordinate space into several regions according to the preset block width and block height, and assigns active grid points located in the same region to the same spatial sub-block; For sub-blocks with fewer than a preset minimum block threshold, they are merged with adjacent normal-sized sub-blocks to ensure that the size of each sub-block meets the solution efficiency requirements of the coherent Ising machine quantum solver.
4. The tokamak plasma current density reconstruction system based on the coherent Ising machine according to claim 3, characterized in that, The partitioning strategy employed by the spatial partitioning processor can be replaced by any of the following strategies: The striping block strategy divides active grid points into blocks of equal length in the order of the global index. The graph partitioning strategy is to divide the graph composed of active grid points into several subgraphs by using the absolute values of the elements of the posterior precision matrix as edge weights, so as to minimize the total weight of the cut edges. An adaptive block partitioning strategy is adopted, which dynamically adjusts the block size based on the number of available bits in the coherent Ising machine quantum solver.
5. The tokamak plasma current density reconstruction system based on the coherent Ising machine according to claim 1, characterized in that, The block coordinate descent iterative controller includes an initialization module, a main update module, a sub-iteration refinement module, a momentum extrapolation module, and a convergence judgment module, wherein: The main update module is used to extract the local objective function matrix and coupling vector of each sub-block, and encode them into a quadratic unconstrained binary optimization sub-problem and submit it to the coherent Ising machine quantum solver for solving. The sub-iterative refinement module is used to fix the integer part of the solution and use the fractional part as a variable to reconstruct and solve the quadratic unconstrained binary optimization sub-problem with higher bit precision based on the current master update solution, so as to improve the accuracy of the solution. The momentum extrapolation module is used to calculate the momentum coefficient based on the momentum strategy of the fast iterative shrinkage threshold algorithm, and to apply extrapolation to the current solution to accelerate convergence.
6. The tokamak plasma current density reconstruction system based on the coherent Ising machine according to claim 1, characterized in that, The block coordinate descent iterative controller also includes a quadratic unconstrained binary optimization encoding module, used to perform the following operations before submitting the subproblem to the coherent Ising machine quantum solver: Estimate the range of values for each continuous variable based on the eigenvalues of the posterior precision matrix; Each continuous variable is encoded using an integer quantum variable, and represented as a weighted sum of multiple binary bits; Substituting the binary code into the local objective function of the sub-block, we can expand it to obtain a quadratic unconstrained binary optimization form. The quadratic unconstrained binary optimization problem is transformed into a standard Ising model through variable transformation, and then solved by the coherent Ising machine quantum solver.
7. The tokamak plasma current density reconstruction system based on the coherent Ising machine according to claim 1, characterized in that, The conditional autoregressive prior model used by the Bayesian posterior modeling processor can be replaced by any of the following models: The advanced squared exponential prior model, which uses a Gaussian process kernel function to describe the spatial correlation of current density; Laplace prior model, used to induce sparsity of solutions; The mixture of Gaussian process prior model, which is composed of multiple Gaussian process priors, is used to capture multi-scale spatial correlations.
8. The tokamak plasma current density reconstruction system based on the coherent Ising machine according to claim 1, characterized in that, The coherent Ising machine quantum solver is a dedicated quantum simulation hardware based on an optical parametric oscillator. When the coherent Ising machine quantum solver is unavailable, the block coordinate descent iterative controller seamlessly switches to a classical solver, which is one or more combinations of fast simulated annealing algorithm, D-Wave quantum annealing machine, branch and bound method or tabu search algorithm.
9. The tokamak plasma current density reconstruction system based on the coherent Ising machine according to claim 1, characterized in that, The external electromagnetic measurement data acquired by the electromagnetic diagnostic data acquisition unit includes signals from magnetic probes, magnetic flux loops, and Roco coils; the system decomposes a large-scale quadratic optimization problem with an overall dimension of hundreds of dimensions into multiple sub-problems with a dimension of tens of dimensions, so that the number of quadratic unconstrained binary optimization bits after encoding each sub-problem does not exceed the upper limit of the hardware bit capacity of the coherent Ising machine quantum solver.