A method of symmetric alternating direction multiplier optimization based on extrapolation acceleration

Abstract

The application provides a symmetrical alternating direction multiplier optimization method based on extrapolation acceleration, which can improve the convergence speed of the algorithm with only a small increase in calculation cost. The extrapolation acceleration mechanism is introduced on the basis of the symmetrical alternating direction multiplier method; firstly, the optimization problem to be solved is decomposed into two sub-problems which can be solved alternately, and an augmented Lagrange function is constructed; secondly, the intermediate multiplier update is introduced after the first original variable update, so that the second original variable can use the updated dual information; then, the extrapolation term is constructed according to the current multiplier update amount and the intermediate multiplier information, and the dual variable is updated at a high speed; finally, by setting the extrapolation step length which satisfies the theoretical upper bound constraint, the numerical stability and convergence of the algorithm are ensured while the convergence speed is improved.

Description

Technical Field

[0001] This invention relates to the field of convex optimization and parameter identification technology, specifically a symmetric alternating direction multiplier optimization method based on extrapolation acceleration. Background Technology

[0002] With the rapid development of large-scale data processing, intelligent control, signal processing, and machine learning, many practical engineering problems can be transformed into convex optimization problems with constraints. For example, in applications such as parameter identification, sparse signal reconstruction, direction-of-arrival estimation, and image denoising, it is often necessary to improve the computational efficiency and convergence speed of optimization algorithms while ensuring solution accuracy. The alternating direction multiplier method is a commonly used constrained optimization solution method. It decomposes complex optimization problems into several easier-to-handle subproblems and alternately updates the original variables and dual multipliers, thereby achieving an efficient solution to large-scale optimization problems. This method has advantages such as clear structure, ease of distributed implementation, and wide applicability, and is therefore widely used in engineering calculations and parameter estimation problems. However, the traditional alternating direction multiplier method has the following problems: (1) The convergence speed is slow in high-precision solutions and large-scale optimization problems, often requiring a large number of iterations, and the computational efficiency needs to be improved.

[0003] (2) The traditional multiplier update method does not make sufficient use of intermediate iteration information. When updating the second original variable, it may rely on older dual information, which affects the overall convergence performance of the algorithm.

[0004] (3) Although the existing symmetric alternating direction multiplier method introduces intermediate multiplier updates, it does not make full use of historical iteration trajectory and extrapolation acceleration information, and there is still room to further improve the convergence speed. Summary of the Invention

[0005] To address the issues of slow convergence speed and numerous iterations in traditional alternating direction multiplier methods for large-scale convex optimization, parameter estimation, sparse signal processing, and image restoration, this invention provides a symmetric alternating direction multiplier optimization method based on extrapolation acceleration, which can improve the algorithm's convergence speed with only a small increase in computational cost.

[0006] The technical solution of this invention is as follows: a symmetric alternating direction multiplier optimization method based on extrapolation acceleration, characterized by comprising the following steps: S1: Construct two separable convex optimization problems to be solved for θ1 and θ2: ; st A1x1+A2x2=b; Where x1 and x2 are variables to be optimized, θ1(x1) and θ2(x2) are convex functions, A1 and A2 are known matrices, and b is a known vector; S2: Based on the separable convex optimization problem, construct the augmented Lagrangian function L: L(x1,x2,λ)=θ1(x1)+θ2(x2)-λ T (A1x1+A2x2-b)+(ρ / 2)·||A1x1+A2x2-b|| 2 ; Where λ is the Lagrange multiplier, ρ>0 is the penalty parameter, and ||·|| 2 The square of the second normal form; S3: Set the initial parameters of the algorithm: Given an initial value x1 0 x2 0 and λ 0 Set the penalty parameter ρ, intermediate multiplier update parameter α, extrapolation step size β, convergence accuracy ε, and maximum number of iterations K. max , where α∈(0,1); S4: In the k-th iteration, fix x2 k and λ k The optimal update value of the first original variable x1 in the current iteration is obtained by minimizing the augmented Lagrange function to update the first variable x1: ; In the formula, L is the augmented Lagrangian function; S5: Based on the updated x1 k+1 For the intermediate multiplier λ k+1 / 2 Update; S6: Fixed x1 k+1 and λ k+1 / 2 Update the second original variable x2; ; S7: Calculate the unextrapolated multiplier update value: based on the updated x1 k+1 and x2 k+1 Calculate the dual variables of the intermediate form ; S8: Utilize the currently unextrapolated multipliers With intermediate multiplier λ 1 / 2+k The difference between them is used to construct extrapolation terms and update the final multipliers; S9: Determine if the algorithm meets the stopping condition; If the stopping condition is met, stop the iteration and proceed to step S10; If the condition is not met, let k = k + 1 and return to step S4 to continue the iteration; S10: Output the optimal solution result; Once the algorithm meets the convergence condition, it outputs the final original variables x1, x2, and the corresponding multiplier variable λ as the solution to the optimization problem.

[0007] Its further features are: The extrapolation step size satisfies the following constraints: 0 < β < (1-α) / α; In step S5, the expression for the intermediate multiplier is: λ k+1 / 2 =λ k -αρ(A1x1 k+1 +A2x2 k -b); In step S7, the expression for the dual variable is: ; In step S8, the expression for the updated final multiplier is: ; Where β is the extrapolation step size; The stopping conditions include: determining whether the algorithm has converged, the change in the objective function being less than a preset threshold, the constraint residual being less than a preset threshold, or the number of iterations reaching the maximum number of iterations K. max ; Where, if ||x i k+1 -x i k If ||<ε, it means the algorithm has converged.

[0008] This application provides a symmetric alternating direction multiplier optimization method based on extrapolation acceleration. It introduces an extrapolation acceleration mechanism on top of the symmetric alternating direction multiplier method. First, the optimization problem to be solved is decomposed into two subproblems that can be solved alternately, and an augmented Lagrangian function is constructed. Second, after updating the first original variable, an intermediate multiplier update is introduced, allowing the second original variable to utilize the updated dual information. Then, an extrapolation term is constructed based on the current multiplier update and the intermediate multiplier information to accelerate the update of the dual variable. Finally, by setting an extrapolation step size that satisfies theoretical bounds, the algorithm maintains numerical stability and convergence while improving convergence speed. This method only needs to sequentially complete the updates of the first original variable, intermediate multipliers, second original variable, and extrapolated multipliers during the iteration process, without complex modifications to the original subproblem structure. Therefore, compared to the traditional alternating direction multiplier method, the additional computational cost is smaller. Furthermore, by introducing an extrapolation acceleration mechanism during the multiplier update process and constraining the extrapolation step size with theoretical bounds, the algorithm effectively reduces the number of iterations and improves convergence speed while ensuring convergence and numerical stability. Therefore, the algorithm used in this invention has good computational efficiency, stability, and application value. Attached Figure Description

[0009] Figure 1 Example diagram of a simulated current signal with an incident angle of 30°; Figure 2 Example diagram of a simulated current signal with an incident angle of 45°; Figure 3 Examples comparing the performance of different algorithms under the same noise conditions; Figure 4 Monte Carlo simulation results using different algorithms; Figure 5 The results of image restoration using different algorithms on a damaged color image; Figure 6 This is a flowchart illustrating the method. Detailed Implementation

[0010] To address the slow convergence speed and numerous iterations of traditional alternating direction multiplier methods in large-scale convex optimization, parameter estimation, sparse signal processing, and image restoration, as well as the insufficient utilization of historical iteration information and intermediate multiplier information in existing symmetric alternating direction multiplier methods, this invention provides an extrapolation-based accelerated symmetric alternating direction multiplier optimization method, namely, Extrapolation-based Accelerated Symmetric Alternating Direction Multiplier Method (EBA-SADMM). This method introduces an extrapolation acceleration mechanism based on the symmetric alternating direction multiplier method, constructing acceleration directions using current and intermediate multiplier update information, thereby improving the algorithm's convergence speed with only a small increase in computational cost. Furthermore, this method provides a theoretical bound on the extrapolation step size, ensuring that the algorithm maintains convergence and numerical stability during accelerated iteration.

[0011] like Figure 6 As shown, this application includes a symmetric alternating direction multiplier optimization method based on extrapolation acceleration, which includes the following steps.

[0012] S1: Construct two separable convex optimization problems to be solved for θ1 and θ2: ; st A1x1+A2x2=b; Where x1∈R n1 and x2∈R n2 Let θ1(x1) and θ2(x2) be the variables to be optimized, and let A1∈R. m×n1 and A2∈R m ×n2 Given a matrix, b∈R m Let n1, n2, and m be known vectors; n1, n2, and m represent the dimensions.

[0013] "st" is an abbreviation for "subject to", meaning "constrained by" or "satisfied with". The goal of the separable convex optimization problem constructed in this application is to find the values ​​of x1 and x2 that minimize θ1(x1) + θ2(x2) under the constraint that x1 and x2 are constrained by the equation A1x1 + A2x2 = b.

[0014] S2: Based on the separable convex optimization problem, construct the augmented Lagrangian function L: L(x1,x2,λ)=θ1(x1)+θ2(x2)-λ T (A1x1+A2x2-b)+(ρ / 2)·||A1x1+A2x2-b|| 2 ; Where, λ∈R mFor Lagrange multipliers, ρ>0 is the penalty parameter, ||·|| 2 It is the square of the second normal form.

[0015] S3: Set the initial parameters of the algorithm: Given an initial value x1 0 x2 0 and λ 0 Set the penalty parameter ρ, intermediate multiplier update parameter α, extrapolation step size β, convergence accuracy ε, and maximum number of iterations K. max , where α∈(0,1). The penalty parameter ρ is used to control the weight of the constraint residuals in the augmented Lagrangian function; the intermediate multiplier update parameter α is used to control the update magnitude of the dual variable.

[0016] The extrapolation step size satisfies the following constraints to ensure algorithm convergence: 0 < β < (1-α) / α.

[0017] S4: In the k-th iteration, fix x2 k and λ k The optimal update value of the first original variable x1 in the current iteration is obtained by minimizing the augmented Lagrange function to update the first variable x1: ; In the formula, L is the augmented Lagrangian function.

[0018] Step S4 is equivalent to performing local optimization on the first variable given the current second variable and multiplier.

[0019] S5: Based on the updated x1 k+1 For the intermediate multiplier λ k+1 / 2 Update; The expression for the intermediate multiplier is: λ k+1 / 2 =λ k -αρ(A1x1 k+1 +A2x2 k -b); The constraint residual information after the first variable update is fed back into the multiplier update in a timely manner, so that the subsequent update of the second variable can utilize the new dual information. Step S5 is used to introduce new dual information between the updates of the two original variables, so that the subsequent variable update can utilize the intermediate multiplier information, thereby improving the problem of lag in multiplier information in traditional ADMM.

[0020] S6: Fixed x1 k+1 and λ k+1 / 2 The second variable x2 is optimized and updated using the information from the intermediate multipliers; .

[0021] In this method, an intermediate multiplier λ is introduced. k+1 / 2 This makes the update of the second original variable x2 no longer completely dependent on the old multiplier information, thereby improving the utilization efficiency of dual information and multiplier information.

[0022] S7: Calculate the updated values ​​of the multipliers without extrapolation; According to the updated x1 k+1 and x2 k+1 Calculate the dual variables of the intermediate form The dual variable represents the temporary multiplier value obtained after updating the two original variables; The expression for the dual variable is: ; Dual variables It reflects the direction of the multiplier correction after the current original variable is updated.

[0023] S8: To further accelerate algorithm convergence, an extrapolation mechanism is introduced to speed up updates: Using the current unextrapolated multipliers With intermediate multiplier λ 1 / 2+k The difference between them is used to construct an extrapolation term, and the final multiplier is updated; the expression for the updated final multiplier is: ; Where β is the extrapolation step size, used to control the acceleration intensity.

[0024] Indicates the current direction of multiplier update; extrapolation term By utilizing the current multiplier update direction, the iteration convergence speed can be improved without significantly increasing computational complexity. Through this extrapolation term, the algorithm can utilize the current and historical iteration trajectory information to further advance along the current iteration trend, thereby reducing the number of iterations required to reach convergence and accelerating convergence while maintaining stability.

[0025] Determine the theoretical range of the extrapolation step size: In order to avoid iterative oscillation or divergence caused by an excessively large extrapolation step size, this method gives the theoretical constraint condition for the extrapolation step size: 0 < β < (1-α) / α; (1-α) / α is the theoretical bound of the extrapolation step size β. When β is too large, although it may enhance the acceleration effect, it may also lead to iteration instability. When β satisfies the above theoretical constraints, the positive definiteness of the correlation matrix can be guaranteed, thereby ensuring that the algorithm has convergence and numerical stability.

[0026] S9: After completing the above update, determine whether the algorithm meets the stopping condition. The specific method is as follows: If any of the following conditions are met, the iteration stops and step S10 is executed; Otherwise, let k = k + 1 and return to step S4 to continue the iteration; Condition 1: Calculate the change in the variable between two consecutive iterations, satisfying ||x i k+1 -x i k ||<ε; Condition 2: The change in the objective function is less than a preset threshold. Condition 3: The constraint residual is less than the preset residual threshold; Condition 4: The number of iterations reaches the maximum number of iterations K. max ; Among them, the preset change threshold, the preset residual threshold, and the maximum number of iterations K are... max Adjust settings according to the actual application scenario.

[0027] S10: Output the optimal solution result; Once the algorithm meets the convergence condition, it outputs the final original variables x1, x2, and the corresponding multiplier variable λ as the solution to the optimization problem.

[0028] Traditional alternating direction multiplier methods (ADMMs) typically obtain the optimal solution by alternately updating the original variables and dual multipliers when solving convex optimization problems with equality constraints. This method is simple in structure and widely applicable, but when facing large-scale parameter estimation, sparse optimization, signal processing, and image restoration problems, it often requires a large number of iterations to achieve high accuracy, resulting in relatively slow convergence. Symmetric alternating direction multiplier methods, by adding intermediate multiplier updates between the updates of two original variables, improve the problem of delayed dual information in traditional ADMMs to some extent, but they still do not fully utilize the historical information and intermediate multiplier information generated during iteration. Therefore, this invention addresses these issues by introducing an extrapolation multiplier update mechanism. This mechanism improves the utilization efficiency of dual information through intermediate multiplier updates and accelerates convergence by using the current iteration direction through extrapolation terms. A theoretical bound on the extrapolation step size is also provided, thereby improving the convergence speed while ensuring the algorithm's convergence and numerical stability. Compared with traditional ADMMs and SADMMs, this invention does not require changes to the basic solution structure of the original subproblems and involves less additional computation. This method can be used for efficient solutions to problems such as large-scale convex optimization, parameter estimation, sparse signal optimization, array signal processing, and image restoration. Theoretical analysis and numerical simulations show that the proposed method has faster convergence speed and better computational efficiency compared to traditional methods in problems such as direction-of-arrival estimation and image denoising.

[0029] Examples of how this method can be applied in actual production work are illustrated below.

[0030] Example 1: LASSO is a class of models widely used in scenarios such as key gene screening in genomics, extraction of core indicators for financial risk control, sparse image denoising and compression, and marketing performance attribution analysis. It is used for automatic feature selection of high-dimensional data, avoiding model overfitting, or constructing interpretable sparse prediction models. LASSO introduces an L1 regularization term to construct an optimization model that combines sparsity and fitting ability. Its objective function is the sum of least squares loss and L1 penalty term. To efficiently solve problems like LASSO, the industry widely adopts ADMM-type methods. These algorithms decompose the original problem into multiple parallel subproblems and alternately update the variables, fully utilizing distributed computing resources and significantly reducing the computational complexity of a single iteration, thereby quickly obtaining high-precision sparse solutions. The key to solving this problem lies in designing a suitable optimization method and efficiently handling the non-smooth subproblems caused by the L1 norm.

[0031] When applying this method to the LASSO parameter estimation problem, for the following sparse optimization model: , st x1-x2=0; In the above formula, A∈R m×n Let b ∈ R be the measurement matrix. m The observation vector is μ>0, which is the regularization parameter. We can set: θ1(x1) = 1 / 2||Ax1-b|| 2 θ2(x2)=μ||x2||1; When applying LASSO to image denoising, b represents the observation vector obtained after vectorizing the noisy image, and A represents the matrix corresponding to the image observation. x1 represents the image vector to be restored obtained from the data fidelity term, and x2 represents the auxiliary image vector used to apply sparse regularization constraints. The equality constraint x1=x2 ensures that the two ultimately correspond to the same restored image.

[0032] The solution is then obtained by substituting the extrapolation acceleration symmetric alternating direction multiplier framework proposed in this invention.

[0033] Example 2: In the algorithm verification stage of radar, sonar, 5G communication and other fields, direction of arrival (DOA) estimation experiments are often carried out. DOA estimation is a core array signal processing technology that uses a sensor array to receive signals and determines the spatial incident angle of the signal source through phase difference or time delay analysis. Experiments typically construct a model including the sensor array, with ADMM and its variant SADMM serving as efficient low-level optimization solvers. In the experiment, incoming wave signals from the 30° and 45° directions are converted into analog current forms that can be collected by the sensor array, providing raw observation inputs that conform to physical characteristics for DOA estimation.

[0034] When this method is applied to the direction-of-arrival (DOA) estimation problem, for DOA estimation models containing sparse noise, it can be transformed into an optimization problem with auxiliary variable constraints: , st vz=0; Where b is the target parameter vector to be estimated, v is the noise vector, and z is the auxiliary variable. By rewriting the problem as two separable structures, this method can be used to solve it, thereby improving the convergence efficiency of parameter estimation. S represents the known source signal matrix, containing prior waveform information. M Let μ represent the M*M identity matrix. μ is the regularization parameter.

[0035] When this method is applied to the direction-of-arrival estimation problem, the specific application process is as follows.

[0036] First, consider the array receiver model containing both Gaussian noise and sparse impulse noise: X = BS + U + V; Where, X∈R M×N Here, M is the array receiving signal matrix, M is the number of array elements, and N is the number of snapshots; X∈R Q×N Given the source signal matrix, containing prior waveform information of the source signals; Q is the number of signals; matrix S is Q*N dimensional; U is the Gaussian noise matrix, and V is the sparse impulse noise matrix; B∈R M×Q Let be the array response matrix to be estimated, and satisfy: B = C(θ)G(g); In the formula, C(θ) is the array manifold matrix composed of the directions of arrival to be estimated, and θ is the direction of arrival vector; G(g) is the source signal complex amplitude matrix, and g is the source signal amplitude vector. Therefore, B contains both direction of arrival and source signal amplitude information. According to the relationship between matrix vectorization and the Kronecker product, we have: ; Therefore, the array receiving model can be transformed into the following vector form: The vectors are defined as follows: x = vec(X), b = vec(B); u=vec(U), v=vec(V); In the formula, vec(·) represents a column vector formed by stacking matrix elements column by column. Represents the Kronecker product; I MLet represent the M×M identity matrix; x is the vectorized result of the array received signal matrix X; b is the vectorized result of the array response matrix B, which is the target parameter vector to be estimated in this model; u and v are the Gaussian noise vector and the sparse noise vector, respectively.

[0037] Assume the vectorized Gaussian noise follows a zero-mean complex Gaussian distribution: ; Where, σ u 2 Let I be the variance of the Gaussian distribution that vector u follows; MN It is an identity matrix with M×N bits; Given b and v, the realized value of Gaussian noise can be obtained from the residuals of the observation model. Therefore, the conditional probability density of the observation vector x can be written as: ; The quadratic term in the conditional probability density is used to measure the deviation between the model-reconstructed signal and the actual received signal, thus constituting the data fidelity term of the optimization model.

[0038] Based on the generalized Gaussian distribution assumption, the prior probability density of sparse noise v can be expressed as: ; Where p(V) is proportional to exp(-k p ||v|| p p ), k p p is the regularization parameter; p controls the tail characteristics and sparsity of the noise distribution. When p=1, the generalized Gaussian distribution degenerates into a Laplace distribution, and its negative log-prior corresponds to the L1 norm regularization term.

[0039] According to Bayes' theorem, the posterior probability of the variable to be estimated satisfies: ; Using the maximum a posteriori estimation criterion, i.e., maximizing the above posterior probability, and taking the negative logarithm of the posterior probability while ignoring constant terms irrelevant to the variable being estimated, we obtain: ; Further multiplying the objective function by a positive constant and unifying the ratio of the two noise parameters as the regularization parameter μ, we can obtain the combined form of the L2 data fidelity term and the p-norm sparse regularization term: ; μ=k p σ u 2 ; Since the distribution parameters of real sparse noise are usually difficult to obtain accurately, μ can be used as an adjustable regularization parameter in practical solutions to balance signal reconstruction accuracy and noise sparsity. For the L1 sparse regularization case, let p=1. To separate the non-smooth regularization term and the quadratic data fidelity term, an auxiliary variable z is introduced, and a consistency constraint z=v is applied, resulting in: ; st vz=0; Here, z is an auxiliary copy of v. The constraint vz=0 ensures that the two remain consistent at the optimal solution, thus not changing the feasible and optimal solutions of the original problem.

[0040] Corresponding to the two optimization models of this invention, b and v form the first variable block, and the auxiliary variable z is used as the second variable block, that is: ; The two corresponding objective functions are defined as follows: ; ; Here, x1 consists of the parameter vector b to be estimated and the sparse noise vector v, and x2 is an auxiliary variable block. The observation vector x = vec(X) is the known array received data, which has a different meaning from the first variable block x1. If the constraint is further written in the standard form A1x1 + A1x2 = c, then we can take: A1=[0 MN×MQ I MN A2 = -I MN c=0; Among them, 0 MN×MQ Represents a 0 matrix of dimension MN×MQ; Thus, the original problem is transformed into an optimization problem in which the objective function is separable and the variables are related by linear equation constraints, satisfying the basic structure of the optimization framework adopted in this invention.

[0041] The hardware platform for the Direction of Arrival (DOA) estimation experiment in this embodiment consists of a computer and an oscilloscope. First, two actual current signals are acquired using the oscilloscope to simulate radar signals, and then imported into the computer to construct a known source signal matrix. The array receiver model is numerically constructed in the computer. In this experiment, the DOAs of the two signal sources are set to 30° and 45°, respectively, and Gaussian noise and sparse impulse noise are added to the observed signals. The experimental software uses MATLAB, implementing the ADMM, SADMM, and EBA-SADMM algorithms respectively. Twenty independent Monte Carlo experiments are conducted under random noise conditions to compare the estimation errors and convergence performance of different algorithms for the target parameter vector and noise vector.

[0042] In the appendix Figures 1-2 middle, Figure 1 and Figure 2 Two analog current signals, respectively, are injected into the sensor from directions of 30° and 45°, and are used to construct observation data in the direction-of-arrival estimation problem. Figure 3 The diagram shows a comparison of error variations under the same noise conditions for the alternating direction multiplier method, the symmetrical alternating direction multiplier method, and the technical solution of this application, namely the symmetrical alternating direction multiplier method based on extrapolation acceleration. Figure 4 The results are Monte Carlo simulations of the alternating direction multiplier method, the symmetric alternating direction multiplier method, and the symmetric alternating direction multiplier method based on extrapolation acceleration in this application.

[0043] Figure 1 and Figure 2 This is the time-domain waveform of two current signals acquired by an oscilloscope. The horizontal axis represents time; the vertical axis represents current amplitude. Each data point on the curve represents the instantaneous current value at the corresponding moment, and the vertical oscillation of the waveform represents the positive or negative change of the current over time. Figure 1 and Figure 2 This represents the input signal used for simulation; the two acquired current signals correspond to different incident angles. The algorithm's parameter estimation capability should refer to... Figure 3 .

[0044] Figure 3 The x-coordinate k represents the number of iterations of the algorithm. The y-coordinate of the upper part is δ. v The lower half of the ordinate δ represents the relative estimation error of the noise vector v. b This represents the relative estimation error of the array parameter vector b to be estimated.

[0045] δ v =||vv ’’ ||2 / ||v ’’ ||2; δ e =||bb ’’ ||2 / ||b ’’ ||2; Where v is the noise vector, v ’’ b is the true value of the noise component; b is the array parameter vector to be estimated. ’’ This represents the actual value of the array parameters; ||·||2 is the L2 normal form.

[0046] Figure 3 The closer the ordinate is to 0, the closer the estimation result is to the true value; the faster the curve decreases, the faster the algorithm converges. The green, blue, and red curves in the figure correspond to ADMM, SADMM, and EBA-SADMM, respectively. The red curve decreases the fastest, indicating that the EBA-SADMM proposed in this application requires fewer iterations under the same error conditions and has a faster convergence speed.

[0047] From the appendix Figure 3 It can be seen that, under the same noise conditions, compared with the traditional alternating direction multiplier method (ADMM) and the symmetric alternating direction multiplier method (SADMM) in the figure, the symmetric alternating direction multiplier method based on extrapolation acceleration proposed in this application (EBA-SADMM) has a smaller estimation error, indicating that this method can maintain good convergence performance and parameter estimation effect even in more complex noise environments.

[0048] based on Figure 1 and Figure 2 The data was used to conduct the Monte Carlo experiment. Figure 4 The x-axis k represents the number of iterations in the algorithm. The y-axis δ in the top row... v The lower row ordinate δ represents the relative estimation error of the noise vector v. b This represents the relative estimation error of the array parameter vector b to be estimated. The three columns correspond to ADMM, SADMM, and EBA-SADMM, respectively. Each thin curve in the figure represents the Monte Carlo simulation result obtained under a single independent random noise condition. Multiple curves represent repeated runs of the algorithm under different random noise samples. The faster the curve decreases, the faster the algorithm converges; the smaller the final error, the higher the estimation accuracy; the more concentrated the different curves are, the less affected the algorithm is by random noise and the better its stability. Based on Figure 4 The data shows that in multiple random noise simulations, the errors of the three algorithms gradually decrease with the increase of the number of iterations. However, EBA-SADMM can achieve a smaller error with fewer iterations, showing more robust convergence characteristics, indicating that it has better convergence speed and resistance to random noise.

[0049] Furthermore, this invention can also be applied to image denoising and image restoration problems. For problems such as robust principal component analysis, low-rank matrix restoration, and sparse noise separation, the low-rank regularization term, sparse regularization term, and constraints in the objective function can be decomposed into different sub-problems, and then the extrapolation of this method is used to accelerate the symmetric alternating direction multiplier method for iterative solution. The hardware platform for the image denoising experiment is a computer used for numerical computation and image processing, without involving external image acquisition equipment. A color image is selected as the original image, and sparse impulse noise is artificially added to the computer to construct a contaminated image. The experimental software uses MATLAB, and an image restoration model is established based on robust principal component analysis, in which the nuclear norm is used to characterize the low-rank structure of the image, and the L1 norm is used to constrain the sparse noise. The model is solved using ADMM, SADMM, and EBA-SADMM respectively, and the restored images and relative errors at the 500th, 1000th, 1500th, and 2000th iterations are recorded to evaluate the image denoising effect and convergence efficiency of different algorithms.

[0050] Figure 5 This is a comparison of the results of image restoration of a damaged color image using three methods: alternating direction multiplier method, symmetric alternating direction multiplier method, and the technical solution of this application. Figure 5 In the text, "Iter=500, 1000, 1500, 2000" represent the image restoration results output when the algorithm runs to the 500th, 1000th, 1500th, and 2000th iterations, respectively. Figure 5 Each row corresponds to an algorithm: the first row is ADMM, the second row is SADMM, and the third row is EBA-SADMM; each column represents the image restoration results of the three algorithms with the same number of iterations. Observing the same row from left to right, we can determine the image restoration process of a certain algorithm; comparing the three rows in the same column, we can determine the restoration performance of the three algorithms with the same number of iterations.

[0051] As the number of iterations increases, the noise in the image gradually decreases, and the outline and details of the person become clearer. At the same number of iterations, if an algorithm produces an image with less noise, clearer outlines, and more complete details, it indicates that the algorithm converges faster and achieves better restoration results. (See attached...) Figure 5 As can be seen, in the image restoration experiment, the symmetric alternating direction multiplier method EBA-SADMM based on extrapolation acceleration proposed in this application can obtain a clearer restored image with fewer iterations, indicating that the extrapolation acceleration mechanism can improve the convergence speed in the image restoration problem.

[0052] In summary, the symmetric alternating direction multiplier optimization method based on extrapolation acceleration proposed in this application can fully utilize intermediate multiplier update information and extrapolation acceleration mechanism to improve the convergence speed and stability of the algorithm. Compared with the traditional alternating direction multiplier method and the symmetric alternating direction multiplier method, the method of this invention has the advantages of smaller error, faster convergence, and better robustness, and can be effectively applied to problems such as parameter estimation, direction of arrival estimation, sparse signal processing, and image restoration.

Claims

1. A symmetric alternating direction multiplier optimization method based on extrapolation acceleration, characterized in that, It includes the following steps: S1: Construct two separable convex optimization problems to be solved for θ1 and θ2: ; st A1x1+A2x2=b; Where x1 and x2 are variables to be optimized, θ1(x1) and θ2(x2) are convex functions, A1 and A2 are known matrices, and b is a known vector; S2: Based on the separable convex optimization problem, construct the augmented Lagrangian function L: L(x1,x2,λ)=θ1(x1)+θ2(x2)-λ T (A1x1+A2x2-b)+(ρ / 2);||A1x1+A2x2-b|| 2 ; Where λ is the Lagrange multiplier, ρ>0 is the penalty parameter, and ||·|| 2 The square of the second normal form; S3: Set the initial parameters of the algorithm: Given an initial value x1 0 x2 0 and λ 0 Set the penalty parameter ρ, intermediate multiplier update parameter α, extrapolation step size β, convergence accuracy ε, and maximum number of iterations K. max , where α∈(0,1); S4: In the k-th iteration, fix x2 k and λ k The optimal update value of the first original variable x1 in the current iteration is obtained by minimizing the augmented Lagrange function to update the first variable x1: ； In the formula, L is the augmented Lagrangian function; S5: Based on the updated x1 k+1 For the intermediate multiplier λ k+1 / 2 Update; S6: Fixed x1 k+1 and λ k+1 / 2 Update the second original variable x2; ； S7: Calculate the unextrapolated multiplier update value: based on the updated x1 k+1 and x2 k+1 Calculate the dual variables of the intermediate form ; S8: Utilize the currently unextrapolated multipliers With intermediate multiplier λ 1 / 2+k The difference between them is used to construct extrapolation terms and update the final multipliers; S9: Determine if the algorithm meets the stopping condition; If the stopping condition is met, stop the iteration and proceed to step S10; If the condition is not met, let k = k + 1 and return to step S4 to continue the iteration; S10: Output the optimal solution result; Once the algorithm meets the convergence condition, it outputs the final original variables x1, x2, and the corresponding multiplier variable λ as the solution to the optimization problem.

2. The symmetric alternating direction multiplier optimization method based on extrapolation acceleration according to claim 1, characterized in that: The extrapolation step size satisfies the following constraints: 0<β< (1-α) / α。 3. The symmetric alternating direction multiplier optimization method based on extrapolation acceleration according to claim 1, characterized in that: In step S5, the expression for the intermediate multiplier is: l k+1 / 2 =λ k -ar(A1 x1 k+1 +A2 x2 k -b)。 4. The symmetric alternating direction multiplier optimization method based on extrapolation acceleration according to claim 1, characterized in that: In step S7, the expression for the dual variable is: 。 5. The symmetric alternating direction multiplier optimization method based on extrapolation acceleration according to claim 1, characterized in that: In step S8, the expression for the updated final multiplier is: ； Where β is the extrapolation step size.

6. The symmetric alternating direction multiplier optimization method based on extrapolation acceleration according to claim 1, characterized in that: The stopping conditions include: determining whether the algorithm has converged, the change in the objective function being less than a preset threshold, the constraint residual being less than a preset threshold, or the number of iterations reaching the maximum number of iterations K. max ; Where, if ||x i k+1 -x i k If ||<ε, it means the algorithm has converged.

Images