Fully automated optimization method for hollow blade section configuration
By introducing a centroid position constraint function into the topology optimization model, the problem of controlling the centroid position of hollow moving blades was solved, realizing fully automatic optimization of hollow moving blades, meeting engineering design requirements and improving structural strength.
Patent Information
- Authority / Receiving Office
- CN · China
- Patent Type
- Patents(China)
- Current Assignee / Owner
- SHANGHAI TURBINE
- Filing Date
- 2021-09-13
- Publication Date
- 2026-06-09
AI Technical Summary
In the existing technology, the topology optimization method for hollow moving blades is difficult to control the centroid position, which makes it difficult to apply the optimization results to actual engineering, and it is difficult to balance the quality and structural strength of hollow moving blades.
By introducing a centroid position constraint function into the topology optimization model, and through finite element analysis and mathematical programming algorithms, the cross-sectional configuration of the hollow moving blade is optimized to ensure a reasonable centroid position, reduce material usage, and guarantee structural strength.
The fully automated optimization of hollow moving blades was achieved, the centroid position met the engineering design requirements, the weight was significantly reduced and the structural strength was improved, overcoming the shortcomings of existing technologies.
Smart Images

Figure CN115809517B_ABST
Abstract
Description
Technical Field
[0001] This invention relates to the field of turbine machinery technology, and in particular to a fully automated method for optimizing the cross-sectional configuration of hollow moving blades. Background Technology
[0002] Turbomachinery is a type of bladed fluid machinery. A rotor equipped with blades rotates at high speed. As fluid flows through the channels between the blades, forces interact between the blades and the fluid, thereby converting energy. The blades, rotating at high speed, bear the enormous centrifugal force generated by their own movement; excessive centrifugal force can cause blade breakage and failure.
[0003] In today's world, where environmental protection, energy conservation, and economic benefits are paramount, turbine users have a strong demand for improving the energy conversion efficiency of turbine machinery. Lengthening the turbine's terminal blades to increase the exhaust area is one effective and feasible solution for improving energy conversion efficiency. Currently, most turbine terminal blades are made of steel. As the blades lengthen, the centrifugal force increases, and the blade length has reached the strength limit of the steel. Reducing the mass of the moving blades while keeping the blade length constant, thereby reducing the centrifugal force and ensuring safe operation, is particularly significant in the design of ultra-long blades.
[0004] Titanium alloys possess advantages such as low density, high specific strength, vibration frequency close to that of steel, and superior resistance to corrosion and water erosion. Using titanium alloys as the material for long, final-stage blades effectively overcomes the blade length limitations imposed by the strength limits of steel. 3D printing is currently one of the commonly used methods for processing titanium alloy parts. As an emerging cutting-edge manufacturing technology, its basic principle is "layer-by-layer manufacturing." Unlike traditional manufacturing methods, 3D printing technology can manufacture parts of arbitrary shapes based on three-dimensional models, without being constrained by the shape of the part and without requiring specialized molds.
[0005] Besides changing the material of the moving blade, adopting a hollow moving blade configuration is also an effective way to reduce the mass of the moving blade. The key to the design of a hollow moving blade structure is the form and size of the cavity and the layout of the internal stiffeners, which have a significant impact on the structural performance of the hollow blade. In existing technologies, the hollow moving blade configuration is generally designed using classical topology optimization methods, which makes it impossible to control the center of gravity of the hollow moving blade, making it difficult to directly apply the optimization results to practical engineering. Summary of the Invention
[0006] In view of the shortcomings of the prior art described above, the technical problem to be solved by the present invention is to provide a fully automatic optimization method for the cross-sectional configuration of hollow moving blades. By adding a centroid position constraint function to the topology optimization model, the centroid position of the optimized configuration of the hollow moving blade meets the design requirements of actual engineering, greatly reducing the mass of the hollow moving blade and ensuring the structural strength of the hollow moving blade.
[0007] To address the aforementioned technical problems, this invention provides a fully automated method for optimizing the cross-sectional configuration of hollow moving blades, comprising the following steps:
[0008] A topology optimization model for the airfoil section of a hollow moving blade is established. The topology optimization model includes design variables, objective function, material usage constraint function, and centroid position constraint function. The design variable is the element pseudo density vector, and the objective function is the structural compliance function.
[0009] The airfoil cross-section is meshed using finite element methods. The finite element elements located at the outer contour edge of the airfoil cross-section are non-design elements, and the pseudo-density ρ of the material in the non-design elements is... j The value is 1, and the remaining finite element elements are design elements. The pseudo-density ρ of the material in the design elements is 1. i For ρ min ≤ρ i ≤1.0.
[0010] After applying load conditions and displacement boundary conditions, and solving the problem using the finite element method, the displacement response value of each node in the finite element mesh is obtained. Furthermore, the sensitivity of at least one of the objective function, material quantity constraint function, and centroid position constraint function to the pseudo-density of the material in the design element is calculated and analyzed.
[0011] A mathematical programming algorithm is selected to solve the topology optimization model, and the optimization solution result of the topology optimization model is obtained;
[0012] Determine whether the optimization solution has converged. If it has, filter the optimization solution to obtain the optimal configuration and perform smoothing design on the optimal configuration. If not, jump back to one of the remaining steps other than the current step until the optimization solution converges.
[0013] Preferably, the unit pseudo-density vector is X = {ρ1, ρ2, ..., ρ...} N} T ,
[0014] Where, ρ i represents the pseudo density of the material in the i-th design unit; N represents the number of design units; T represents the inverse and transpose, or the transpose and inverse.
[0015] Preferably, the objective function is F = Compliance = U T KU,
[0016] Where U represents the displacement vector; K represents the overall structural stiffness matrix; and T represents the inverse and transpose, or the transpose and inverse.
[0017] Preferably, the centroid position constraint function is:
[0018] Where, r i ={x i ,y i} T Let r0 = {x0, y0} represent the centroid coordinates of the i-th design unit. T The target centroid coordinates of the airfoil section are represented, N represents the number of design elements, and ρ represents the number of design elements. i Let ζ represent the pseudo density of the material in the i-th design unit, and ζ be the centroid constraint parameter.
[0019] Preferably, the material usage constraint function is:
[0020] Where, ρ i V represents the pseudo-density of the material in the i-th design unit. i ρ represents the area of the i-th design unit. i V i V represents the material usage of the i-th design unit, V0 is the total area of all design units, Δ represents the volumetric constraint parameter of the material usage, and N represents the number of design units.
[0021] Preferably, the filtering process uses the Heaviside formula.
[0022] The Heaviside formula is
[0023] in, ρ represents the filtered and corrected topological variables. i Let represent the pseudo density of the material in the i-th design unit, β be the correction factor, and e be the natural constant, which is usually taken as e≈2.718.
[0024] Preferably, the sensitivity analysis formula for the objective function is as follows:
[0025] in, This represents the sensitivity of the objective function to the pseudo-density of the material in the design element. This represents the derivative of the overall structural stiffness matrix with respect to the material pseudo-density of the design element. This represents the derivative of the load vector with respect to the material pseudo-density of the design element.
[0026] Preferably, the fully automatic optimization method for the cross-sectional configuration of the hollow moving blade further includes: providing the outer contour of the airfoil cross-section, the material elastic modulus, Poisson's ratio, load conditions, and displacement boundary conditions of the hollow moving blade.
[0027] Preferably, the step of establishing a topology optimization model for the airfoil cross-section of the hollow moving blade, wherein the topology optimization model includes design variables, objective function, material usage constraint function, and centroid position constraint function, the design variable being the element pseudo-density vector, and the objective function being the structural compliance function, includes: setting the volume fraction constraint parameter of the material usage, and obtaining the target centroid coordinates of the airfoil cross-section.
[0028] Preferably, the mathematical programming algorithm is one of sequential linear programming, sequential quadratic programming, and moving asymptote method.
[0029] As described above, the fully automated optimization method for the cross-sectional configuration of hollow moving blades of the present invention has the following beneficial effects: Based on topology optimization theory, the present invention provides a fully automated optimization method for the cross-sectional configuration of hollow moving blades to achieve an innovative fully automated design of hollow moving blade configurations for turbine machinery. The entire process is automated without human intervention, and the resulting optimal configuration has a reasonable centroid position. This avoids the drawback of designers being unable to control the centroid position when manually adjusting the topology of the moving blade cavity, thus reducing material usage while ensuring the structural strength of the moving blade. The fully automated optimization method of the present invention introduces a design variable of material pseudo-density into a given design domain. Unlike existing topology optimization methods, the fully automated optimization method of the present invention adds a centroid position constraint function for the airfoil cross-section of the hollow moving blade to the topology optimization model to fully consider the design requirements of the centroid position of the hollow moving blade. This ensures that the final airfoil cross-section configuration has a reasonable centroid position, overcoming the drawback of existing topology optimization methods where the centroid position of the hollow moving blade cannot be controlled, making it difficult to directly apply the hollow moving blade to practical engineering. Therefore, the fully automatic optimization method for the cross-sectional configuration of the hollow moving blade of the present invention, by adding a centroid position constraint function to the topology optimization model, makes the centroid position of the optimized configuration of the hollow moving blade meet the design requirements of actual engineering, greatly reducing the mass of the hollow moving blade and ensuring the structural strength of the hollow moving blade. Attached Figure Description
[0030] Figure 1 A schematic diagram showing the target centroid position of the airfoil section of a hollow moving blade;
[0031] Figure 2 This is a schematic diagram showing the setting of the material pseudo-density for the design unit;
[0032] Figure 3 The diagram shows the design domain and non-design domain of the airfoil section;
[0033] Figure 4 The diagram shows a variable density penalty function.
[0034] Figure 5 The diagram shows the filtering function.
[0035] Figure 6 A schematic diagram showing the finite element mesh model of the airfoil section;
[0036] Figure 7 This is a schematic diagram showing the optimal configuration of the airfoil section;
[0037] Figure 8 The diagram shows the final configuration after smoothing the optimal airfoil section design.
[0038] Figure 9 The flowchart shown is an embodiment of the fully automatic optimization method for the cross-sectional configuration of hollow moving blades according to the present invention. Detailed Implementation
[0039] The following specific embodiments illustrate the implementation of the present invention. Those skilled in the art can easily understand other advantages and effects of the present invention from the content disclosed in this specification.
[0040] It should be understood that the structures, proportions, sizes, etc., depicted in the accompanying drawings of this specification are merely for illustrative purposes to aid those skilled in the art and are not intended to limit the scope of the invention. Therefore, they have no substantial technical significance. Any modifications to the structure, changes in proportions, or adjustments to size, without affecting the effectiveness and purpose of the invention, should still fall within the scope of the technical content disclosed in this invention. Furthermore, the terms such as "upper," "lower," "left," "right," "middle," and "one" used in this specification are merely for clarity and are not intended to limit the scope of the invention. Changes or adjustments to their relative relationships, without substantially altering the technical content, should also be considered within the scope of the invention's implementation.
[0041] like Figure 1 , Figure 2 as well as Figure 3 As shown, this invention provides a fully automated optimization method for the cross-sectional configuration of hollow moving blades, comprising the following steps:
[0042] A topology optimization model for the airfoil section of a hollow moving blade is established. The topology optimization model includes design variables, objective function, material usage constraint function, and centroid position constraint function. The design variable is the element pseudo density vector, and the objective function is the structural compliance function.
[0043] The aforementioned airfoil cross-section is meshed using finite element methods. The finite element elements located at the outer contour edge of the airfoil cross-section are non-design elements, and the pseudo-density ρ of the material in the non-design elements is... j The value is 1, and the remaining finite element elements are design elements. The pseudo-density ρ of the material in the design elements is 1. i For ρmin ≤ρ i ≤1.0.
[0044] After applying load conditions and displacement boundary conditions, and solving the problem using the finite element method, the displacement response value of each node in the finite element mesh is obtained. Furthermore, the sensitivity of at least one of the objective function, material quantity constraint function, and centroid position constraint function to the pseudo-density of the material in the design element is calculated and analyzed.
[0045] A mathematical programming algorithm is selected to solve the above topology optimization model, and the optimized solution result of the topology optimization model is obtained.
[0046] Determine whether the optimization solution has converged. If it has, filter the optimization solution to obtain the optimal configuration and perform smoothing design on the optimal configuration. If not, jump back to one of the remaining steps other than the current step until the optimization solution converges.
[0047] Based on topology optimization theory, this invention provides a fully automated optimization method for the cross-sectional configuration of hollow moving blades, enabling an innovative configuration design for hollow moving blades in turbine machinery. This fully automated optimization method introduces a design variable representing the pseudo-density of the material within a given design domain. Unlike existing topology optimization methods, this invention's fully automated optimization method, to fully consider the design requirements for the centroid position of the hollow moving blade (when the hollow moving blade is made of homogeneous material, its centroid coincides with its center of mass), adds a centroid position constraint function for the airfoil cross-section of the hollow moving blade to the topology optimization model. This ensures that the final airfoil cross-section configuration has a reasonable centroid position, overcoming the drawback of existing topology optimization methods where the centroid position of the hollow moving blade cannot be controlled, making it difficult to directly apply the hollow moving blade to practical engineering.
[0048] Specifically, the fully automated optimization method of this invention includes the following steps: establishing a topology optimization model of the airfoil cross-section of the hollow moving blade, the topology optimization model including design variables, objective function, material usage constraint function, and centroid position constraint function; performing finite element mesh generation on the above airfoil cross-section (the number of meshes is controlled between 2000 and 10000 to provide the node coordinates, finite element numbers, and center coordinates of the finite element elements for subsequent finite element analysis), typically the resulting finite element is a quadrilateral isoparametric element. Based on the variable density method, the pseudo-density of the material in the quadrilateral isoparametric element constitutes the elements of the above design variables. See also... Figure 3 The outer contour edge of the blade section is defined as the non-design domain, and the area of the blade section outside the outer contour edge is defined as the design domain. The finite element elements located at the outer contour edge of the blade section are non-design elements. Considering the functional requirements such as the shape of the hollow moving blade, the non-design elements are located within the non-design domain. The pseudo-density ρ of the material of the non-design elements is...j =1, ρ j =1.0(j=J1,J2,L,J) M ), where J represents the row number of the non-design element and M represents the number of non-design elements; the remaining finite element elements are design elements, located within the design domain, and the pseudo-density ρ of the material in the design element. i For ρ min ≤ρ i ≤1.0. See also Figure 2 ρ = 1 indicates that the design unit is filled with solid material, and ρ = 0 indicates that the design unit is empty and is considered void. To avoid singularities in the overall structural stiffness matrix, ρ min A value of 0.001 can be used. Apply load conditions and displacement boundary conditions, and after solving using the finite element method, obtain the displacement response value of each node in the finite element mesh. Calculate and analyze the sensitivity of at least one of the objective function, material usage constraint function, and centroid position constraint function to the pseudo-density of the material in the design element. Select a mathematical programming algorithm to solve the above topology optimization model and obtain the optimization solution result. Determine whether the optimization solution result has converged. If it has, filter the optimization solution result to obtain the optimal configuration, and perform smoothing design on the obtained optimal configuration to obtain the final configuration of the blade section. If not, jump back to one of the remaining steps except the current step until the optimization solution result converges.
[0049] Therefore, the fully automatic optimization method for the cross-sectional configuration of the hollow moving blade of the present invention, by adding a centroid position constraint function to the topology optimization model, makes the centroid position of the optimized configuration of the hollow moving blade meet the design requirements of actual engineering, greatly reducing the mass of the hollow moving blade and ensuring the structural strength of the hollow moving blade.
[0050] The pseudo-density vector of the above unit is X = {ρ1, ρ2, L, ρ N} T ,
[0051] Where, ρ i represents the pseudo density of the material in the i-th design unit; N represents the number of design units; T represents the inverse and transpose, or the transpose and inverse.
[0052] The objective function above is F = Compliance = U T KU,
[0053] Where U represents the displacement vector; K represents the overall structural stiffness matrix; and T represents the inverse and transpose, or the transpose and inverse.
[0054] The centroid position constraint function is as follows:
[0055] Where, ri ={x i ,y i} T Represents the centroid coordinates of the i-th design unit; see [link / reference] Figure 1 r0 = {x0, y0} T The target centroid coordinates of the airfoil section are represented, N represents the number of design elements, and ρ represents the number of design elements. i Let ξ represent the pseudo density of the material in the i-th design unit, and ζ be the centroid constraint parameter, typically ξ = 0.0001 to 0.01.
[0056] The above material usage constraint function is:
[0057] Where, ρ i V represents the pseudo-density of the material in the i-th design unit. i ρ represents the area of the i-th design unit. i V i V represents the material usage of the i-th design unit, V0 is the total area of all design units, Δ represents the volume ratio constraint parameter of the material usage, which is generally 0.3 to 0.6, and N represents the number of design units.
[0058] The above filtering process uses the Heaviside formula (i.e., the filtering function).
[0059] The Heaviside formula is
[0060] in, ρ represents the filtered and corrected topological variables. i Let represent the pseudo-density of the material in the i-th design unit, β be the correction factor (usually 50), and e represent the natural constant (usually e ≈ 2.718). When β is 50, the filtering function is as follows: Figure 5 As shown.
[0061] The sensitivity analysis formula for the above objective function is as follows:
[0062] in, This represents the sensitivity of the objective function to the pseudo-density of the material in the design element. This represents the derivative of the overall structural stiffness matrix with respect to the material pseudo-density of the design element. This represents the derivative of the load vector with respect to the material pseudo-density of the design element. Therefore, the fully automatic optimization method of this invention can perform topology optimization design for different working conditions, such as pressure loads, inertial loads, and gravity loads, and designers can apply load conditions according to the actual working conditions.
[0063] The fully automated optimization method for the aforementioned hollow moving blade cross-sectional configuration further includes providing the blade profile, material elastic modulus, Poisson's ratio, load condition P, and displacement boundary conditions of the hollow moving blade. Furthermore, based on the concept of variable density, the fully automated optimization method of this invention can also use a variable density penalty function to link the pseudo-density of the material in the design unit with the material elastic modulus of the hollow moving blade. The variable density penalty function adopts a power function form, i.e. Where E i Let represent the elastic modulus of the i-th design unit, E represent the material elastic modulus of the hollow moving blade, and η represent the penalty factor. When η is 3, 4, or 5, the variable density penalty function is as follows: Figure 4 As shown.
[0064] The above-mentioned steps for establishing a topology optimization model for the airfoil section of a hollow moving blade include design variables, objective function, material usage constraint function, and centroid position constraint function. The design variable is the element pseudo density vector, and the objective function is the structural compliance function. Specifically, these steps include: setting the volume fraction constraint parameter of the material usage, and obtaining the target centroid coordinates of the airfoil section.
[0065] The mathematical programming algorithm described above is one of the following: sequential linear programming, sequential quadratic programming, and moving asymptote method.
[0066] like Figure 6 and Figure 7 As shown, the airfoil cross-section of a hollow moving blade is used as an optimization case: the optimization objective is selected as compliance under static load, and the material usage within the design domain does not exceed 35% of the entire airfoil cross-section. The entire airfoil cross-section is divided into 4409 four-node isoparametric elements. Among them, the number of non-design elements located at the outer contour is 872. The moving asymptotic method is used to solve this topology optimization model, and the optimal configuration obtained is shown in [reference]. Figure 7 See the final configuration after smoothing design. Figure 8 It can be seen that the optimal configuration of the obtained blade section is a truss-like structure, which is one of the configurations with the highest stiffness. The force transmission path is reasonable, and the expected design requirements are met.
[0067] In addition, see Figure 9 As an embodiment of the fully automated optimization method of the present invention, the method includes the following steps:
[0068] S1, Establish the initial topology optimization model of the airfoil section;
[0069] S2, set the material volume ratio, and obtain the target centroid coordinates of the blade section;
[0070] S3 uses four-node isoparametric elements for mesh generation;
[0071] S4, set the pseudo-density of materials for non-design units and design units respectively;
[0072] S5, Set the centroid constraint parameter ξ;
[0073] S6, apply the material's elastic modulus, Poisson's ratio, etc.;
[0074] S7, apply load and displacement boundary conditions;
[0075] S8, Finite element method, calculates the sensitivity of the objective function and constraint function for each design element;
[0076] S9, Select an algorithm to solve;
[0077] S10: Determine whether the optimization solution has converged. If yes, proceed to S11; otherwise, jump back to S6.
[0078] S11, perform smoothing design on the optimal configuration;
[0079] S12, End.
[0080] In summary, the fully automated optimization method for the hollow moving blade cross-section configuration of the present invention, by adding a centroid position constraint function to the topology optimization model, ensures that the centroid position of the optimized configuration of the hollow moving blade meets the design requirements of actual engineering, greatly reducing the mass of the hollow moving blade and ensuring its structural strength. Therefore, the present invention effectively overcomes the various shortcomings of the prior art and has high industrial application value.
[0081] The above embodiments are merely illustrative of the principles and effects of the present invention and are not intended to limit the invention. Any person skilled in the art can modify or alter the above embodiments without departing from the spirit and scope of the present invention. Therefore, all equivalent modifications or alterations made by those skilled in the art without departing from the spirit and technical concept disclosed in the present invention should still be covered by the claims of the present invention.
Claims
1. A fully automated optimization method for the cross-sectional configuration of a hollow moving blade, characterized in that, Includes the following steps: A topology optimization model for the airfoil cross-section of a hollow moving blade is established. The topology optimization model includes design variables, an objective function, material usage constraint functions, and centroid position constraint functions. The design variables are element pseudo-density vectors, and the objective function is a structural compliance function. The element pseudo-density vector is... ;in, The pseudo density of the material in the i-th design unit is represented by N; the number of design units is represented by T; and the inverse is represented by the transpose after inversion, or vice versa. The objective function is: ,in, Represents the displacement vector; represents the overall stiffness matrix of the structure; T represents the transpose after inversion, or the inversion after transposition; The sensitivity analysis formula for the objective function is as follows: ,in, This represents the sensitivity of the objective function to the pseudo-density of the material in the design element. This represents the derivative of the overall structural stiffness matrix with respect to the material pseudo-density of the design element. This represents the derivative of the load vector with respect to the material pseudo-density of the design element. Represents the displacement vector; The centroid position constraint function is: ,in, This represents the centroid coordinates of the i-th design unit. The coordinates of the target centroid of the airfoil section are given, and N represents the number of design elements. This represents the pseudo-density of the material in the i-th design unit. These are the centroid constraint parameters; The airfoil cross-section is meshed using finite element methods. The finite element elements located at the outer contour edge of the airfoil cross-section are non-design elements, and the pseudo-density ρ of the material in the non-design elements is... j The value is 1, and the remaining finite element elements are design elements. The pseudo-density ρ of the material in the design elements is 1. i for ; After applying load conditions and displacement boundary conditions, and solving the problem using the finite element method, the displacement response value of each node in the finite element mesh is obtained. Furthermore, the sensitivity of at least one of the objective function, material quantity constraint function, and centroid position constraint function to the pseudo-density of the material in the design element is calculated and analyzed. A mathematical programming algorithm is selected to solve the topology optimization model, and the optimization solution result of the topology optimization model is obtained; Determine whether the optimization solution has converged. If it has, filter the optimization solution to obtain the optimal configuration and perform smoothing design on the optimal configuration. If not, jump back to one of the remaining steps other than the current step until the optimization solution converges.
2. The fully automatic optimization method for the cross-sectional configuration of hollow moving blades according to claim 1, characterized in that: The material usage constraint function is: , in, This represents the pseudo-density of the material in the i-th design unit. This represents the area of the i-th design unit. This represents the material usage of the i-th design unit. The total area of all design units. The volumetric constraint parameter represents the material usage, and N represents the number of design units.
3. The fully automatic optimization method for the cross-sectional configuration of hollow moving blades according to claim 1, characterized in that: The filtering process uses the Heaviside formula, which is: , in, These are the topological variables after filtering correction. This represents the pseudo-density of the material in the i-th design unit. The correction factor is e, which represents the natural constant and is usually taken as e ≈ 2.
718.
4. The fully automatic optimization method for the cross-sectional configuration of hollow moving blades according to claim 1, characterized in that: The fully automatic optimization method for the cross-sectional configuration of the hollow moving blade also includes: providing the outer contour of the airfoil cross-section, the material elastic modulus, Poisson's ratio, load conditions, and displacement boundary conditions of the hollow moving blade.
5. The fully automatic optimization method for the cross-sectional configuration of hollow moving blades according to claim 1, characterized in that: The steps of establishing a topology optimization model for the airfoil section of the hollow moving blade, which includes design variables, objective function, material usage constraint function, and centroid position constraint function, with the design variable being the element pseudo-density vector and the objective function being the structural compliance function, include: setting the volume fraction constraint parameter of the material usage and obtaining the target centroid coordinates of the airfoil section.
6. The fully automatic optimization method for the cross-sectional configuration of hollow moving blades according to claim 1, characterized in that: The mathematical programming algorithm is one of the following: sequential linear programming, sequential quadratic programming, and moving asymptote method.