Expensive constrained multi-objective optimization method and system based on constrained hyper-volume expectation

Abstract

The invention discloses an expensive constrained multi-objective optimization method and system based on constrained hyper-volume expectation. The method comprises the following steps: firstly, initializing a population and training an agent model; then calculating a correlation coefficient between the solution with the best population convergence and the constraint violation degree of the solution, and performing classification processing; then, optimization is carried out on the agent model, and a new solution is selected by utilizing a filling criterion based on constraint hyper-volume expectation to carry out evaluation and is filled into a file; further selecting and screening individuals by using an improved environment to enter a next-generation population; and finally, judging whether a stop condition is reached or not, if so, outputting a final population, otherwise, turning to execute an iteration process. According to the method, the optimization efficiency is remarkably improved, local optimum is avoided, constraint conditions are precisely processed, and the population can be more effectively converged to a better feasible region.

Description

Technical Field

[0001] The present invention belongs to the technical field of multi-objective optimization in computer technology, and relates to an adaptive expensive-constrained multi-objective optimization method and system, and particularly to an adaptive expensive-constrained multi-objective optimization method and system based on the improvement of the expected constrained hypervolume, which can be applied to engineering design, scientific computing, and physical simulation. Background Art

[0002] In the fields of engineering design, scientific computing, physical simulation, etc., expensive-constrained multi-objective optimization problems often occur. Such problems refer to a class of optimization problems where the computational costs of the objective functions and constraint conditions are very high. Such problems have the following characteristics:

[0003] (1) Multi-objective: There are multiple objective functions, and trade-offs need to be made among them. Usually, there is no clear single optimal solution, but a set of non-dominated solutions (Pareto solution set), which makes the problem more complex.

[0004] (2) Constraint conditions: The optimization problem also needs to satisfy certain constraint conditions.

[0005] (3) High cost: The evaluation costs of the objective functions and constraints are very high, and expensive simulations or experiments may be required. Therefore, a large number of evaluations cannot be carried out like traditional optimization problems.

[0006] Due to the high evaluation cost, a key challenge in solving expensive-constrained multi-objective optimization problems is how to ensure the effectiveness of the optimization process and effectively handle the objectives and constraint conditions within a limited computational budget. Therefore, sample-based optimization methods, such as efficient global optimization algorithms, are usually required to approximate the optimization process and reduce the need for evaluating the objective functions and constraint functions. Currently, one of the main methods for solving expensive-constrained multi-objective optimization problems is the surrogate model-assisted evolutionary algorithm, which reduces the computational overhead by approximating the surrogate model to replace the expensive objective functions and constraint functions, and combines the evolutionary algorithm to efficiently explore the search space. However, the existing expensive-constrained multi-objective algorithms still have the following problems:

[0007] (1) The constraint surrogate model may mislead the algorithm search because the surrogate model cannot accurately capture the true form of complex non-linear constraints, especially near the constraint boundary, and may produce large prediction errors, resulting in the optimization process deviating from the true feasible region.

[0008] (2) The algorithm is prone to falling into local optima. During the optimization process, when the surrogate model guides the search direction, it focuses more on the information within the current local region and lacks continuous and effective exploration of the global search space. When the algorithm enters a local relatively optimal region, it is difficult to jump out and misses the region where the global optimal solution is located, greatly limiting the performance of the algorithm.

[0009] (3) Lack of an effective model filling strategy. The reason is that existing methods often fail to reasonably select the data points for updating the surrogate model, resulting in uneven data distribution, insufficient filling or overfitting in some regions, thus affecting the global search ability of the optimization process, wasting computing resources, and possibly leading to an imbalance between exploration and exploitation.

[0010] Therefore, how to solve the problem of the surrogate model misleading the algorithm, help the algorithm jump out of the local optimal region, and reasonably select the data points for updating the surrogate model is a key technical problem that urgently needs to be solved. Summary of the Invention

[0011] To solve the above technical problems of expensive constraints and multi-objectives in the fields of engineering design, scientific computing, physical simulation, etc., the present invention provides an adaptive expensive constraint multi-objective optimization method and system based on the improvement of constrained hypervolume expectation.

[0012] The technical solution adopted by the method of the present invention is: an expensive constraint multi-objective optimization method based on constrained hypervolume expectation, including the following steps:

[0013] Step 1: Initialize the population and train the surrogate model;

[0014] Step 2: Calculate the correlation coefficient between the solution with the best population convergence and its constraint violation degree, and classify and process it;

[0015] Step 3: Perform optimization on the surrogate model;

[0016] Step 4: Use the filling criterion based on constrained hypervolume expectation to select new solutions for evaluation and fill them into the archive;

[0017] Step 5: Use improved environmental selection to screen individuals into the next generation population;

[0018] Step 6: Determine whether the stop condition is reached. If so, output the final population; otherwise, return to execute Steps 2 to 5.

[0019] Preferably, in Step 1, the Latin hypercube sampling method is used to generate the initial population; the specific implementation includes the following sub-steps:

[0020] Step 1.1: Define the decision variables x1, x2,..., x

[0022] ,

[0021] , i , , i , , d and their value ranges [a i , b i , i = 1,..., d;

[0021] Step 1.2: Divide the value range of each variable into N sub-intervals with equal probability, where N is the population size;

[0022] Step 1.3: Independently and randomly select points within each variable sub-interval, and determine the point values through \(x\) i,j = a i +(j - 1)Δx i + r i,j Δx i where \(r\) i,j is a uniformly distributed random number in [0, 1], and Δx i is the interval length; \(a\) i represents the lower limit of the value range of the \(i\)-th decision variable \(x\) i , and \(j\) represents the serial number of dividing the sub-interval in Step 1.2;

[0023] Step 1.4: Combine the selected points of each variable to form \(N\) \(d\)-dimensional individuals, constituting the initial population Population.

[0024] Preferably, in Step 1, after sampling, the initial population is evaluated using the true objective function and constraint function; the sampled decision values, objective function, and decision function are used to train the surrogate model.

[0025] The surrogate model regards the unknown function as an implementation of a Gaussian process \(y(x)\), and its mathematical model is:

[0026] \(y(x)=f(x)+\epsilon\);

[0027] wherein, is Gaussian noise with a uniform distribution, and \(\mu\) and \(\sigma\) are constants independent of \(x\);

[0028] Given the collected data \(X = \{x_1, x_2, \ldots, x\) N}\} and its corresponding outputs \(Y = \{y(x_1), y(x_2), \ldots, y(x\) N )\}, for the untested point \(z\) * , the predicted mean * \(\hat{y}(z\) and the predicted variance \(\sigma\) 2 [\hat{f}(z\) * )] are estimated according to the following formula:

[0029]

[0030]

[0031]

[0032] where \(\psi(X)\) is the mean vector of \(X\), and \(k\) * is the covariance between \(X\) and \(z\) *The covariance vector, where \(K\) is the covariance matrix of \(X\); \(k(x, x')\) is the covariance function that measures the similarity between two data points \(x\) and \(x'\), \(I\) represents the identity matrix, and \(\theta\), \(l\) are two scale parameters used to control the covariance magnitude and decay rate. is the noise variance;

[0033] The hyperparameters are determined by maximizing the log marginal likelihood function:

[0034]

[0035] where \(m(X)\) is the predicted mean vector on the training data \(X\), \(N\) is the number of training data points, and \(K\) is the covariance matrix. is the noise variance.

[0036] Preferably, in step 2, the metric \(I\) in the IBEA algorithm ∈+ is used to calculate the solution with the best population convergence; specifically, the convergence metric \(C(x)\) is calculated as follows:

[0037]

[0038]

[0039] where \(x\), \(y\) are any two different solutions in the population, \(P\) is the population, \(g(x)\) is the objective function, and \(m\) represents the number of objective functions;

[0040] The convergence metric \(C(x)\) of the population is calculated, and the top \(X\) solutions are selected in descending order as the solutions with the best population convergence; where \(X\) is a preset value;

[0041] Next, the constraint violation degree \(CV(x)\) of these \(X\) solutions is calculated, and its calculation method is as follows:

[0042]

[0043] where \(g\) i (x) is the \(i\)-th constraint function with respect to \(x\);

[0044] Next, the Pearson correlation coefficient \(r\) between \(C(x)\) and \(CV(x)\) is calculated to measure the correlation between the algorithm convergence and the constraint violation degree; its calculation method is as follows:

[0045]

[0046] where is the mean of the convergence metric \(C(x)\) of the selected \(X\) solutions, represents the mean of the constraint violation degree \(CV(x)\) of the selected \(X\) solutions, and \(n\) is the number of solutions; \(C(x i ), \(CV(x i) represents the convergence degree and constraint violation degree of the \(i\)-th solution \(x\) among the \(X\) solutions taken out i ;

[0047] The value range of the Pearson correlation coefficient \(r\) is \([-1, 1]\); if \(r > 0\) and its absolute value is greater than the threshold \(A\), the changes of \(C(x)\) and \(CV(x)\) are strongly positively correlated, then unconstrained optimization is adopted; if \(r < 0\) and its absolute value is greater than the threshold \(B\), the changes of \(C(x)\) and \(CV(x)\) are strongly negatively correlated, then constrained optimization is adopted; otherwise, multi-objective optimization is carried out by taking the constraint violation degree \(CV(x)\) as an additional objective to achieve a trade-off.

[0048] Preferably, in step 2, when the algorithm is in the early stage, that is, when the number of evaluations is less than the threshold \(C\), the unconstrained optimization strategy is fixedly used to enable the algorithm to explore the search space before using the correlation coefficient; in the later stage of the algorithm, that is, when the number of evaluations is greater than the threshold \(D\), the constrained optimization strategy is fixedly used.

[0049] Preferably, the specific implementation of step 3 includes the following sub-steps:

[0050] Step 3.1: Select the current population as the initial population for optimization;

[0051] Step 3.2: Use the surrogate model to evaluate the objective function value and constraint violation degree of each individual in the initial population; calculate the fitness value of the individual according to the objective function value and constraint violation degree;

[0052] Step 3.3: Perform non-dominated sorting on the individuals in the population; when there is no constraint, directly perform non-dominated sorting according to the objective function value; when considering constraints, first divide the feasible and infeasible solutions, the feasible solutions are sorted by the conventional non-dominated sorting, the infeasible solutions are sorted by the constraint violation degree, and the feasible solutions are given priority; in the case of weak correlation, the \(CV(x)\) value is used as an additional objective for non-dominated sorting; divide the individuals into different non-dominated levels; the individuals within the same level are non-dominated with each other, and there is a domination relationship between different levels;

[0053] Step 3.4: Perform crossover and mutation operations on the selected parent individuals to generate offspring individuals; the crossover operation exchanges the genes of the parent individuals, and the mutation operation randomly perturbs the genes of the individuals to increase the diversity of the population;

[0054] Step 3.5: Combine the parent population and the offspring population into a new population

[0055] Step 3.6: Perform non-dominated sorting and reference point association on the combined population, and select a certain number of individuals to form the next generation population according to the occupancy of the reference point and the non-dominated level of the individuals;

[0056] Step 3.7: Determine whether the termination condition is satisfied. If the termination condition is satisfied, the algorithm ends; otherwise, return to Step 3.2 to continue the iteration.

[0057] Preferably, the specific implementation of Step 4 includes the following sub-steps:

[0058] Step 4.1: After obtaining the optimization result on the surrogate model, cluster the objective values of the algorithm; use the K-Means algorithm to divide the population into μ clusters, where μ is the number of solutions sampled each time;

[0059] Step 4.2: For each cluster, calculate the CEHVI value CEHVI(x) of the solutions therein;

[0060] The specific calculation process is as follows:

[0061] (1) Sample and matrix preparation; When the algorithm is initialized, use Latin hypercube sampling to generate N samples in the region [lb, ub] to form a sample set S, and at the same time create an N × N matrix T to store the dominance relationship between samples; where lb and ub are preset values; M Generate N E samples to form a sample set S, and at the same time create an N E ×N E matrix T to store the dominance relationship between samples; where lb and ub are preset values;

[0062] (2) Normalization processing; Before calculating EHVI, perform a normalization operation on the objectives of the current population and the current Pareto front PF; for the j-th objective, use a simple min-max normalization method to normalize the mean and variance calculated by the surrogate model, and the calculation is as follows:

[0063]

[0064] where, represents the predicted value of the decision variable x on the j-th objective function; represents the maximum value of the j-th objective function that has been evaluated in the archive; represents the minimum value of the j-th objective function that has been evaluated in the archive; represents the predicted variance of the decision variable x on the j-th objective function;

[0065] (3) Update matrix T; Update matrix T according to the current Pareto front PF, remove the rows and columns corresponding to the samples dominated by the current Pareto front PF, and based on the updated matrix T nd , obtain HVI(q i ) by counting rows or columns; where q i is one of the samples;

[0066]

[0067] Among them, is the hypervolume of the set PF;

[0068] (4) Calculate EHVI in different cases; if the variance of an individual is zero, it means that the individual has been evaluated by the expensive function, and at this time EHVI is equal to zero; otherwise, based on N nd uniform samples, calculate EHVI according to the following formula:

[0069]

[0070] where q ij represents the j-th objective of sample q i ; N nd is the number of samples from the uniform distribution and not dominated by the current Pareto front PF; M is the number of objective functions; respectively represent the standardized predicted value and predicted variance of solution x on the j-th objective; φ is the cumulative distribution function;

[0071] (5) Calculate the CEHVI value CEHVI(x) of the solutions therein;

[0072] CEHVI(x) = EHVI(x) × PoF(x);

[0073] Calculate the probability PoF(x) of the solution satisfying the constraint through the degree of constraint violation calculated by the surrogate model;

[0074]

[0075] where, is the j-th predicted constraint value of solution x, and are its predicted constraint mean and uncertainty respectively, P(·) and Φ(·) represent probability and Gaussian cumulative distribution function respectively; p represents the number of constraint functions;

[0076] Step 4.3: From each cluster, select the solution with the largest CEHVI value and add it to the candidate solution set.

[0077] Step 4.4: Select μ candidate solutions for evaluation and put them into the evaluated solution set, and update the model using the evaluated solution set.

[0078] Preferably, the specific implementation of step 5 includes the following sub-steps:

[0079] Step 5.1: Calculate the constraint violation value of each solution, and divide all solutions into two categories: feasible solutions and infeasible solutions according to whether the value is 0 or greater than 0;

[0080] Step 5.2: Extract the objective values of the feasible solutions, perform non-dominated sorting, obtain the dominance levels and the maximum dominance level of each feasible solution, and screen out the non-dominated feasible solutions;

[0081] Step 5.3: Extract the objective values and constraint violation values of the infeasible solutions, traverse each infeasible solution, determine whether it is dominated by the non-dominated feasible solutions, and mark the infeasible solutions that are not dominated;

[0082] Step 5.4: Perform non-dominated sorting on the dominated feasible solutions;

[0083] Step 5.5: Perform non-dominated sorting on the remaining infeasible solutions;

[0084] Step 5.6: Sort the population into a new population in the order of non-dominated feasible solutions, infeasible solutions with better convergence, dominated feasible solutions, and remaining infeasible solutions;

[0085] Step 5.7: Retain a specified number of solutions.

[0086] The technical solution adopted by the system of the present invention is: an expensive constrained multi-objective optimization system based on constrained hypervolume expectation, including:

[0087] One or more processors;

[0088] A storage device for storing one or more programs, which, when executed by the one or more processors, cause the one or more processors to implement the expensive constrained multi-objective optimization method based on constrained hypervolume expectation.

[0089] Compared with the prior art, the beneficial effects of the present invention include:

[0090] (1) The present invention calculates the Pearson correlation coefficient r between the population convergence index C and the constraint violation degree CV, and adaptively selects the optimization strategy according to the positive or negative and absolute value size of r, which is beneficial to saving computing resources.

[0091] (2) The present invention calculates the probability of feasibility PoF of the solution and multiplies it by the quickly calculated EHVI to obtain the constraint-aware CEHVI. The K-Means algorithm is used to cluster the objective values of the NSGA-III optimization results, and the solution with the largest CEHVI value is selected from each cluster and added to the candidate solution set for evaluating and updating the model to prevent the algorithm from falling into local optimality.

[0092] (3) Aiming at the problem that the traditional non-dominated sorting combined with the constraint handling technology overemphasizes the feasibility of solutions, which may mislead the population convergence in the problem of discontinuous feasible regions, the present invention improves the environmental selection strategy. It reduces the priority of the dominated feasible solutions and guides the population to converge to a better feasible region. Brief Description of the Drawings

[0093] The following uses examples and specific implementation manners to further illustrate the technical solution of the present invention. In addition, some drawings are also used in the process of describing the technical solution. For those skilled in the art, without creative efforts, other drawings and the intention of the present invention can also be obtained based on these drawings.

[0094] Figure 1 It is the flowchart of the method for this embodiment. Specific implementation manner

[0095] To facilitate the understanding and implementation of the present invention by those of ordinary skill in the art, the present invention will be further described in detail below in conjunction with the drawings and embodiments. It should be understood that the embodiments described herein are only used to illustrate and explain the present invention, and are not used to limit the present invention.

[0096] Please see Figure 1 , a method for expensive constrained multi-objective optimization based on constrained hypervolume expectation provided in this embodiment includes the following steps:

[0097] Step 1: Initialize the population and train the surrogate model;

[0098] In one implementation manner, the Latin hypercube sampling method is used to generate the initial population. The specific steps are as follows:

[0099] 1) Identify the decision variables x1, x2, …, x d and their value ranges [a i , b i (i = 1, …, d).

[0100] In practical problems, x1, x2, …, x d respectively represent the decision variables after the problem to be solved is modeled as an expensive constrained multi-objective optimization problem.

[0101] In the field of engineering design, taking an engineering design example of an automotive component as an example, this example needs to optimize the design of an energy absorption component to find the optimal combination of the thicknesses of nine walls to improve its impact resistance. The goal of the design optimization is to maximize energy absorption and minimize the compression force, while satisfying the material quality and being small. The decision variables are the thicknesses of nine walls (x1, x2, …, x9), and the range is [0.2, 3], with the unit of millimeters.

[0102] 2) Divide the value range of each variable into N (population size) sub-intervals with equal probability.

[0103] 3) Independently and randomly select points within each variable sub-interval, and through x i,j = a i +(j - 1)Δx i + r i,j Δxi (r i,j is a uniform random number in [0,1], and Δx i is the interval length) to determine the point value. a i represents the lower limit of the value range of the i-th decision variable x i , and j represents the serial number of dividing the sub-interval in step 2);

[0104] 4) Combine the selected points of each variable to form N d-dimensional individuals to constitute the initial population Populaiton.

[0105] After sampling, use the true objective function and constraint function to evaluate the initial population; use the sampled decision values, objective function, and decision function to train a Kriging model.

[0106] The Kringing model regards the unknown function as an implementation of the Gaussian process y(x), and its mathematical model is:

[0107] y(x) = f(x) + ∈

[0108] In the formula, is Gaussian noise with a uniform distribution, and μ and σ are constants independent of x.

[0109] Given the collected data X = {x1, x2,..., x N} and its corresponding output Y = {y(x1), y(x2),..., y(x N ), for the untested point z * , the mean value * and variance σ and variance σ 2 of the approximate value y(z

[0110]

[0111] are estimated according to the following formula: * where ψ(X) is the mean vector of X, and k * is the covariance vector between X and z

[0112]

[0113] I represents the identity matrix, and θ and l are two scale parameters used to control the covariance size and decay rate, is the noise variance;

[0114] The predicted mean is directly used as the predicted value of y(z * ), and the predicted variance σ 2 [f(z* ) reflects uncertainty. The hyperparameters in the formula are determined by maximizing the log marginal likelihood function:

[0115]

[0116] where m(X) is the predicted mean vector on the training data X, N is the number of training data points, K is the covariance matrix, is the noise variance.

[0117] Step 2: Calculate the correlation coefficient between the solution with the best population convergence and its degree of constraint violation, and classify and process it;

[0118] In one implementation, for the current population, select the solution with the best population convergence. The specific selection method uses the metric I in the IBEA algorithm ∈+ to calculate. The specific convergence metric C(x) is calculated as follows:

[0119]

[0120]

[0121] where x and y are any two different solutions in the population, P is the population, f(x) is the objective function, and m represents the number of objective functions; is a binary metric. In the convergence metric C(x), use the summation symbol to calculate the I metric operation between all solutions y in the population P except x and x, and then accumulate to finally obtain the value of C(x).

[0122] Calculate the convergence degree C(x) of the population, and select the top 20 solutions as the solutions with the best population convergence by descending order. Next, calculate the degree of constraint violation CV(x) of these twenty solutions, and its calculation method is as follows:

[0123]

[0124] where g i (x) is the i-th constraint function with respect to x.

[0125] Next, calculate the Pearson correlation coefficient between C(x) and CV(x) to measure the correlation between the algorithm convergence and the degree of constraint violation. Its calculation method is as follows:

[0126]

[0127] where, the mean of the convergence metrics C(x) of the X solutions taken, represents the mean of the degrees of constraint violation CV(x) of the X solutions taken, and n is the number of solutions; C(x i)、CV(x i ) represents the convergence degree and constraint violation degree of the \(i\)-th solution \(x\) among the \(X\) solutions taken out i ;

[0128] The value range of the Pearson correlation coefficient \(r\) is \([-1, 1]\). If \(r > 0\), it indicates that there is a positive linear correlation between \(C(x)\) and \(CV(x)\), and the convergence direction and the direction of constraint optimization may be the same. If \(r < 0\), it indicates that there is a negative linear correlation between \(C(x)\) and \(CV(x)\), and the convergence direction and the direction of constraint optimization may be opposite. The larger the absolute value of \(r\), the stronger the correlation.

[0129] To sum up, when \(r > 0\) and the absolute value is relatively large, here we take \(r > 0.6\), the changes of \(C(x)\) and \(CV(x)\) show a strong positive correlation. In order to save computing resources, unconstrained optimization can be directly adopted. When \(r < 0\) and the absolute value is relatively large, here we take \(r < 0.2\), the changes of \(C(x)\) and \(CV(x)\) show a strong negative correlation, which means that the changes of convergence and constraint violation degree are opposite, and constrained optimization needs to be adopted. When the absolute value of \(r\) is small, the changes of \(C(x)\) and \(CV(x)\) show a weak correlation, which indicates that there is no significant correlation between the two changes. At this time, in order to balance the constraints and the convergence of the objective value, the constraint violation degree \(CV(x)\) is used as an additional objective for multi-objective optimization to achieve a trade-off.

[0130] Specifically, when the algorithm is in the early stage, that is, when the number of evaluations is small, since the algorithm has not fully converged and the exploration of the search space is limited, the correlation coefficient may mislead the selection of the algorithm strategy. Therefore, in the early stage of the algorithm, the unconstrained optimization strategy is fixed to make the algorithm explore the search space and then use the correlation coefficient. In the later stage of the algorithm, the algorithm has fully explored the search space and needs to conduct a fine search at the constraint boundary. Therefore, in the late stage of the algorithm, that is, when the remaining number of evaluations is small, the optimization strategy considering constraints is fixed.

[0131] Step 3: Optimize on the surrogate model;

[0132] In one implementation, NSGA-III is used as an optimizer to optimize on the surrogate model; NSGA-III is a powerful evolutionary algorithm for solving multi-objective optimization problems. Compared with traditional multi-objective optimization algorithms, NSGA-III has better convergence and distribution. When dealing with expensive constraint multi-objective optimization problems, due to the high computational cost of actually evaluating the objective function and constraint conditions, in this embodiment, a surrogate model is used to approximate the real objective function and constraint conditions. NSGA-III can efficiently search for a set of non-dominated solutions with good distribution and convergence on the surrogate model, providing rich choices for subsequent decision-making.

[0133] The specific optimization process of NSGA-III is as follows:

[0134] 1) Initialize the population: Here, the current population is selected as the initial population for optimization.

[0135] 2) Fitness evaluation: Use the surrogate model to evaluate the objective function value and constraint violation degree of each individual in the initial population. Calculate the fitness value of the individual according to the objective function value and constraint violation degree.

[0136] 3) Non-dominated sorting: Perform non-dominated sorting on the individuals in the population. When there are no constraints, directly perform non-dominated sorting according to the objective function value. When considering constraints, first divide the feasible and infeasible solutions. The feasible solutions are sorted by the conventional non-dominated sorting, and the infeasible solutions are sorted by the constraint violation degree, with feasible solutions having priority. In the case of weak correlation, the CV value is used as an additional objective for non-dominated sorting. Divide the individuals into different non-dominated levels. The individuals within the same level are non-dominated with each other, and there is a dominance relationship between different levels.

[0137] 4) Crossover and mutation operations: Perform crossover and mutation operations on the selected parent individuals to generate offspring individuals. The crossover operation exchanges the genes of the parent individuals, and the mutation operation randomly perturbs the genes of the individuals to increase the diversity of the population.

[0138] 5) Merge the populations: Merge the parent population and the offspring population into a new population.

[0139] 6) Environmental selection: Perform non-dominated sorting and reference point association on the merged population. Select a certain number of individuals to form the next generation population according to the occupancy of the reference point and the non-dominated level of the individuals.

[0140] 7) Termination condition judgment: Judge whether the termination condition is satisfied, such as reaching the maximum number of iterations, convergence of the objective function value, etc. If the termination condition is satisfied, the algorithm ends; otherwise, return to step 2 to continue the iteration.

[0141] Step 4: Select new solutions for evaluation using the filling criterion based on the expected constrained hypervolume and fill them into the archive;

[0142] In one implementation, the adopted filling criterion comprehensively considers the constraint satisfaction probability of the solution and the hypervolume expected improvement EHVI, and uses importance sampling to simplify the calculation and combines the constraints to obtain the EHVI considering the constraints. The specific calculation process is as follows.

[0143] In multi-objective optimization, the calculation method of the hypervolume improvement is calculated as follows:

[0144]

[0145] Among them, is the hypervolume of the set PF. The calculation of EHVI can be written as:

[0146]

[0147] where V nd is the non-dominated solution set. However, EHVI has no analytical formula that can be directly used for calculation, and Monte Carlo sampling is usually used, resulting in low computational efficiency. By reconstructing the integral formula with importance sampling, it can be rewritten, and finally, using uniform sampling, the formula for EHVI(x) can be rewritten as:

[0148]

[0149] where N nd is the number of samples from the uniform distribution that are not dominated by the current Pareto front (PF), and q i is one of the samples.

[0150] The specific calculation process is as follows:

[0151] 1) Sample and matrix preparation: When the algorithm is initialized, Latin Hypercube Sampling (LHS) is used to generate N M samples in the region [lb, ub] E to form a sample set S. At the same time, an N E ×N E matrix T is created to store the dominance relationship between samples. lb is usually slightly less than 0, and ub is slightly greater than 1 to ensure that uniform sampling can cover all regions not dominated by PF.

[0152] 2) Normalization: Before calculating EHVI, normalize the objectives of the current population and the current PF. For the j-th objective, use a simple min-max normalization method to normalize the mean and variance calculated by the surrogate model, and the calculation is as follows:

[0153]

[0154] where represents the predicted value of the decision variable w on the j-th objective function; represents the maximum value of the j-th objective function that has been evaluated in the archive; represents the minimum value of the j-th objective function that has been evaluated in the archive; represents the predicted variance of the decision variable x on the j-th objective function;

[0155] 3) Update matrix T: Update matrix T according to the current PF, remove the rows and columns corresponding to the samples dominated by the current PF, and based on the updated matrix T nd , obtain HVI(q i)。

[0156] 4) Calculate EHVI in different cases. If the variance of an individual is zero, it means that the individual has been evaluated by the expensive function, and at this time EHVI is equal to zero; otherwise, based on N nd uniform samples, calculate EHVI according to the following formula:

[0157]

[0158] where q ij represents the j-th objective of sample q i ; N nd is the number of samples from the uniform distribution that are not dominated by the current Pareto front PF; M is the number of objective functions; respectively represent the standardized predicted value and predicted variance of solution x on the j-th objective; φ is the cumulative distribution function;

[0159] The above is the calculation process of EHVI. To consider constraints, we calculate the probability of solution satisfying constraints PoF by considering the degree of constraint violation calculated by the surrogate model. The specific calculation method is as follows:

[0160]

[0161] where is the j-th predicted constraint value of solution x, and are its predicted constraint mean and uncertainty respectively, P(·) and Φ(·) represent probability and Gaussian cumulative distribution function respectively; p represents the number of constraint functions;

[0162] In summary, the calculation method of EHVI considering constraints is obtained:

[0163] CEHVI(x) = CEHVI(x) × PoF(x).

[0164] The filling strategy of the algorithm is specifically implemented as follows:

[0165] 1) After obtaining the optimization result of NSGA-III on the surrogate model, in order to prevent the algorithm from falling into local optimum, first cluster the objective values of the algorithm. Use the K-Means algorithm to divide the population into μ clusters, where μ is the number of solutions sampled each time.

[0166] 2) For each cluster, calculate the CEHVI value of the solutions in it.

[0167] 3) From each cluster, select the solution with the largest CEHVI value and add it to the candidate solution set.

[0168] Through this filling strategy, the algorithm selects μ candidate solutions for evaluation and puts them into the evaluated solution set. Update the model using the evaluated solution set.

[0169] Step 5: Use the improved environmental selection to screen individuals into the next generation population;

[0170] When the traditional non-dominated sorting is combined with the constraint handling technique, the principle of feasibility first is often adopted. That is, it is considered that feasibility is prior. If CV(x1) < CV(x2), then x1 dominates x2. However, this environmental selection strategy overemphasizes the feasibility of the solution, which will cause the feasible solutions to be retained in the population for a long time. However, in the optimization problem with discontinuous feasible regions, retaining some dominated feasible solutions for a long time may mislead the population and prevent it from converging to a better feasible region. Therefore, in one implementation, an improved environmental selection strategy is adopted to reduce the priority of the dominated feasible solutions. The specific process is as follows:

[0171] 1) Divide feasible and infeasible solutions: Calculate the constraint violation value of each solution. According to whether the value is 0 or greater than 0, all solutions are divided into two categories: feasible solutions and infeasible solutions.

[0172] 2) Find non-dominated feasible solutions: Extract the objective values of the feasible solutions, perform non-dominated sorting, obtain the domination levels and the maximum domination level of each feasible solution, and screen out the non-dominated feasible solutions.

[0173] 3) Find infeasible solutions with better convergence: Extract the objective values and constraint violation values of the infeasible solutions, traverse each infeasible solution, judge whether it is dominated by the non-dominated feasible solutions, and mark the infeasible solutions that are not dominated.

[0174] 4) For the dominated feasible solutions, perform non-dominated sorting.

[0175] 5) For the remaining infeasible solutions, perform non-dominated sorting.

[0176] 6) Sort the population into a new population in the order of non-dominated feasible solutions, infeasible solutions with better convergence, dominated feasible solutions, and remaining infeasible solutions.

[0177] 7) Retain a specified number of solutions.

[0178] Step 6: Judge whether the stop condition is reached. If so, output the final population; otherwise, return to execute Steps 2 to 5.

[0179] This embodiment also provides an expensive constraint multi-objective optimization system based on the constrained hypervolume expectation, including:

[0180] One or more processors;

[0181] A storage device for storing one or more programs, which when executed by the one or more processors, cause the one or more processors to implement the expensive constraint multi-objective optimization method based on the constrained hypervolume expectation.

[0182] The present invention will be further described through specific experiments below.

[0183] As shown in Table 1, two current advanced algorithms were selected in this experiment for comparison on the constrained multi-objective optimization MW test set. The problem dimension was 10, and the limited number of evaluations was 300 times. The selected test index was IGD+. Thirty repeated experiments were carried out on the test set to eliminate randomness, and the results given in the table were presented in the form of mean (variance). In addition, the rank-sum test results were also included in the table. In the given results, "+" may indicate significantly better than, "-" indicates significantly worse than, and "=" indicates no significant difference.

[0184] IGD + (Inverted Generational Distance Plus) index is an important measure for evaluating the quality of the non-dominated solution set. This index is improved on the basis of the traditional IGD index. By introducing a modified distance measure, it can more effectively reflect the approximation degree between the non-dominated solution set and the reference point set. It is calculated as follows:

[0185]

[0186] where a i(k) is the nearest objective vector of the reference point z k in the objective vector set A, |Z| is the cardinality of the reference point set Z, and d + is the modified distance, and the calculation method is:

[0187] It can be seen from Table 1 that the method proposed by the present invention performs optimally in most problems of the test set. And through statistical tests, the algorithm proposed by the present invention has significant advantages in most problems.

[0188] Table 1

[0189]

[0190]

[0191] The present invention has the following advantages:

[0192] 1. The optimization efficiency is significantly improved: Through the Kriging surrogate model, the direct evaluation times of expensive objective functions and constraint functions are greatly reduced. In the face of complex expensive optimization problems, compared with traditional algorithms, the present invention can perform effective searches with less computational budget. For example, in the constrained multi-objective traditional test set MW, in the case of 10 dimensions, the present invention can obtain better optimization results with 300 evaluations.

[0193] 2. Avoiding being trapped in local optima: The combination of the K-Means algorithm and the constraint-based EHVI filling strategy can effectively prevent the algorithm from being trapped in local optima. K-Means clusters the optimization results of NSGA-III, guiding global exploration, balancing search and exploitation, and adapting to the diversity of solutions, enabling the algorithm to select optimal solutions from different regions. The constraint-based EHVI comprehensively considers objective optimization and constraint satisfaction, quantifies the potential of solutions, and handles uncertainties, preferentially selecting solutions that are likely to improve the objective and satisfy the constraints. Their collaboration provides a clear orientation for the algorithm in complex multi-objective and constraint environments, enhancing the global search ability and reducing the risk of being trapped in local optima.

[0194] 3. Precisely handling constraint conditions: By calculating the correlation coefficient to adaptively select the optimization strategy, the constraint conditions can be reasonably handled at different stages. An unconstrained optimization strategy is adopted in the early stage of the algorithm to fully explore the search space, and in the later stage, it focuses on fine-grained search of the constraint boundaries. At the same time, the improved environmental selection strategy reduces the priority of dominated feasible solutions, enabling the population to converge more effectively to a better feasible region.

[0195] It should be understood that the above-described embodiments are part of the embodiments of the present invention, rather than all embodiments. Additionally, the technical features in each embodiment or individual embodiment provided by the present invention can be combined with each other arbitrarily to form a feasible technical solution. This combination is not restricted by the order of steps and / or the structural composition mode, but must be based on what can be achieved by those of ordinary skill in the art. When the combination of technical solutions results in contradictions or cannot be implemented, it should be considered that such a combination of technical solutions does not exist and is not within the scope of protection required by the present invention.

[0196] It should be understood that the above description of the preferred embodiments is relatively detailed, and thus it should not be considered as a limitation on the protection scope of the patent of the present invention. Those of ordinary skill in the art, under the inspiration of the present invention and without departing from the scope protected by the claims of the present invention, can still make substitutions or modifications, which all fall within the protection scope of the present invention. The scope of protection requested by the present invention shall be subject to the appended claims.

Claims

1. An expensive constrained multi-objective optimization method based on the constrained hypervolume expectation, characterized in that, It includes the following steps: Step 1: Initialize the population and train the surrogate model; Step 2: Calculate the correlation coefficient between the solution with the best population convergence and its constraint violation degree, and classify and process it; Step 3: Perform optimization on the surrogate model; Step 4: Use the filling criterion based on the constrained hypervolume expectation to select new solutions for evaluation and fill them into the archive; Step 5: Use the improved environmental selection to screen individuals into the next generation population; Step 6: Judge whether the stopping condition is reached. If so, output the final population; otherwise, loop back to execute Steps 2 to 5.

2. The expensive constrained multi-objective optimization method based on the constrained hypervolume expectation according to claim 1, wherein: In Step 1, the Latin hypercube sampling method is used to generate the initial population; the specific implementation includes the following sub-steps: Step 1.1: Define the decision variables x1, x2, …, x for the optimization problem d and their value ranges [a i , b i , where i = 1, …, d; Step 1.2: Divide the value range of each variable into N sub-intervals with equal probability, where N is the population size; Step 1.3: Independently and randomly select points within each variable sub-interval, and determine the point value through x i,j = a i +(j - 1)Δx i + r i,j Δx i where r i,j is a uniform random number in [0, 1], and Δx i is the interval length; a i represents the lower limit of the value range of the i-th decision variable x i and j represents the serial number of the sub-interval division in Step 1.2; Step 1.4: Combine the selected points of each variable to form N d-dimensional individuals, constituting the initial population Populaiton.

3. The expensive constrained multi-objective optimization method based on the constrained hypervolume expectation according to claim 2, characterized in that: In Step 1, after sampling, use the true objective function and constraint function to evaluate the initial population; use the sampled decision values, objective function, and decision function to train the surrogate model; The surrogate model regards the unknown function as an implementation of the Gaussian process y(w), and its mathematical model is: y(x) = f(w) + ∈; where is Gaussian noise with a uniform distribution, and μ and σ are constants independent of x; Given the collected data \(X = \{x_1, x_2, \ldots, x\) N \}\) and its corresponding outputs \(Y=\{y(x_1), y(x_2), \ldots, y(x\) N )\}, for an untested point \(z\) * , the predicted mean * \(\hat{y}(z\) ) and the predicted variance \(\sigma\) 2 ^2[\hat{f}(z\) * )] are estimated by the following formulas: where ψ(X) is the mean vector of X, k * is the covariance vector between X and z * , K is the covariance matrix of X; k(x, x′) is the covariance function that measures the similarity between two data points x and x′, I represents the identity matrix, θ and l are two scale parameters used to control the covariance magnitude and decay rate, is the noise variance; The hyperparameters are determined by maximizing the log marginal likelihood function: where \(m(X)\) is the predicted mean vector on the training data \(X\), \(N\) is the number of training data points, \(K\) is the covariance matrix, and \(\sigma^2\) is the noise variance.

4. The expensive constrained multi-objective optimization method based on the constrained hypervolume expectation according to claim 1, characterized in that: In step 2, the index I in the IBEA algorithm is used ∈+ to calculate the solution with the best population convergence; specifically, the convergence index C(x) is calculated as follows: where x and y are any two different solutions in the population, P is the population, f(x) is the objective function, and m represents the number of objective functions; Calculate the convergence index C(x) of the population, and select the top X solutions as the solutions with the best population convergence by descending order; where X is a preset value; Next, calculate the constraint violation degree CV(x) of these X solutions, and its calculation method is as follows: where g i (x) is the i-th constraint function with respect to x, and m is the number of constraint functions; Next, calculate the Pearson correlation coefficient r between C(x) and CV(x) to measure the correlation between the algorithm convergence and the constraint violation degree; its calculation method is as follows: Among them, The mean value of the convergence index C(x) of the X solutions taken, represents the mean value of the constraint violation degree CV(x) of the X solutions taken, where n is the number of solutions; C(x i ), CV(x i ) represent the convergence degree and the constraint violation degree of the i-th solution x i among the X solutions taken; The value range of the Pearson correlation coefficient r is [-1, 1]; if r > 0 and the absolute value is greater than the threshold A, the changes of C(x) and CV(x) show a strong positive correlation, then unconstrained optimization is adopted; if r < 0 and the absolute value is greater than the threshold B, the changes of C(x) and CV(x) show a strong negative correlation, and constrained optimization is adopted; otherwise, use the constraint violation degree CV(x) as an additional objective for multi-objective optimization to seek a trade-off.

5. The expensive constrained multi-objective optimization method based on the constrained hypervolume expectation according to claim 1, characterized in that: In Step 2, when the algorithm is in the early stage, that is, the number of evaluations is less than the threshold C, the unconstrained optimization strategy is fixed to make the algorithm explore the search space before using the correlation coefficient; in the later stage of the algorithm, that is, the number of evaluations is greater than the threshold D, the constrained optimization strategy is fixed.

6. The expensive constrained multi-objective optimization method based on the constrained hypervolume expectation according to claim 1, wherein The specific implementation of Step 3 includes the following sub-steps: Step 3.1: Select the current population as the initial population for optimization; Step 3.2: Use the surrogate model to evaluate the objective function value and constraint violation degree of each individual in the initial population; calculate the fitness value of the individual according to the objective function value and constraint violation degree; Step 3.3: Perform non-dominated sorting on the individuals in the population; when there are no constraints, directly perform non-dominated sorting according to the objective function values; when considering constraints, first divide the solutions into feasible and infeasible solutions. For feasible solutions, perform conventional non-dominated sorting, and for infeasible solutions, sort them according to the degree of constraint violation, with feasible solutions taking precedence; in the case of weak correlation, use the CV(x) value as an additional objective for non-dominated sorting; divide the individuals into different non-dominated levels; individuals within the same level are non-dominated with respect to each other, and there is a dominance relationship between different levels; Step 3.4: Perform crossover and mutation operations on the selected parent individuals to generate offspring individuals; the crossover operation exchanges the genes of the parent individuals, and the mutation operation randomly perturbs the genes of the individuals to increase the diversity of the population; Step 3.5: Combine the parent population and the offspring population into a new population; Step 3.6: Perform non-dominated sorting and reference point association on the combined population, and select a certain number of individuals to form the next generation population according to the occupancy of the reference points and the non-dominated levels of the individuals; Step 3.7: Determine whether the termination condition is satisfied. If the termination condition is satisfied, the algorithm ends; otherwise, return to Step 3.2 to continue the iteration.

7. The expensive constrained multi-objective optimization method based on the constrained hypervolume expectation according to claim 1, wherein The specific implementation of Step 4 includes the following sub-steps: Step 4.1: After obtaining the optimization results on the surrogate model, cluster the objective values of the algorithm; use the K-Means algorithm to divide the population into μ clusters, where μ is the number of solutions sampled each time; Step 4.2: For each cluster, calculate the CEHVI value CEHVI(x) of the solutions therein; The specific calculation process is as follows: (1) Sample and matrix preparation; when the algorithm is initialized, Latin hypercube sampling is used in the region [lb, ub] M to generate N E samples, forming a sample set S. At the same time, an N E ×N E matrix T is created to store the dominance relationship between samples; where lb and ub are preset values; (2) Normalization processing; before calculating the EHVI, perform normalization operations on the objectives of the current population and the current Pareto front PF; for the j-th objective, use the simple min-max normalization method to normalize the mean and variance calculated by the surrogate model, and the calculation is as follows: Among them, represents the predicted value of the decision variable x on the j-th objective function; represents the maximum value of the j-th objective function evaluated in the archive; represents the minimum value of the j-th objective function evaluated in the archive; represents the predicted variance of the decision variable x on the j-th objective function; (3) Update matrix T; update matrix T according to the current Pareto front PF, remove the rows and columns corresponding to the samples dominated by the current Pareto front PF, and based on the updated matrix T hd , obtain the HVI(q i ) by counting the rows or columns; where q i is one of the samples; Among them, is the hypervolume of the set PF; (4) Calculate EHVI in different cases; if the variance of an individual is zero, it means that the individual has been evaluated by the expensive function, and at this time EHVI is equal to zero; otherwise, calculate EHVI according to the following formula based on N nd uniform samples: where q ij represents the j-th objective of sample q i ; N nd is the number of samples from the uniform distribution that are not dominated by the current Pareto front PF; M is the number of objective functions; respectively represent the standardized predicted value and the prediction variance of solution x on the j-th objective; φ is the cumulative distribution function; (5) Calculate the CEHVI value CEHVI(x) of the solutions therein; CEHVI(x) = EHVI(x) × PoF(x); Calculate the probability PoF(x) of the solution satisfying the constraint through the degree of constraint violation calculated by the surrogate model; where is the j-th predicted constraint value of the solution x, and are its predicted constraint mean and uncertainty respectively, P(·) and Φ(·) represent probability and Gaussian cumulative distribution function respectively; p represents the number of constraint functions; Step 4.3: From each cluster, select the solution with the largest CEHVI value and add it to the candidate solution set; Step 4.4: Select μ candidate solutions for evaluation and put them into the evaluated solution set, and update the model using the evaluated solution set.

8. The expensive constrained multi-objective optimization method based on the constrained hypervolume expectation according to any one of claims 1-7, characterized in that The specific implementation of Step 5 includes the following sub-steps: Step 5.1: Calculate the constraint violation value of each solution, and divide all solutions into two categories: feasible solutions and infeasible solutions according to whether the value is 0 or greater than 0; Step 5.2: Extract the objective values of the feasible solutions, perform non-dominated sorting, obtain the dominance levels and the maximum dominance level of each feasible solution, and screen out the non-dominated feasible solutions; Step 5.3: Extract the objective values and constraint violation values of the infeasible solutions, traverse each infeasible solution, and determine whether it is dominated by the non-dominated feasible solutions, and mark the non-dominated infeasible solutions; Step 5.4: For the dominated feasible solutions, perform non-dominated sorting; Step 5.5: For the remaining infeasible solutions, perform non-dominated sorting; Step 5.6: Sort the population into a new population in the order of non-dominated feasible solutions, infeasible solutions with better convergence, dominated feasible solutions, and remaining infeasible solutions; Step 5.7: Retain a specified number of solutions.

9. An expensive constrained multi-objective optimization system based on the expected value of the constrained hypervolume, characterized in that Comprising: One or more processors; A storage device for storing one or more programs, which when executed by the one or more processors, cause the one or more processors to implement the expensive constrained multi-objective optimization method based on constrained hypervolume expectation according to any one of claims 1 to 8.

Images