CN109814391B - Singular optimal control simultaneous solving method based on partial moving finite element nodes - Google Patents
Singular optimal control simultaneous solving method based on partial moving finite element nodes Download PDFInfo
- Publication number
- CN109814391B CN109814391B CN201910120094.9A CN201910120094A CN109814391B CN 109814391 B CN109814391 B CN 109814391B CN 201910120094 A CN201910120094 A CN 201910120094A CN 109814391 B CN109814391 B CN 109814391B
- Authority
- CN
- China
- Prior art keywords
- finite element
- optimal control
- singular
- node
- control curve
- Prior art date
- Legal status (The legal status is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the status listed.)
- Active
Links
Images
Landscapes
- Management, Administration, Business Operations System, And Electronic Commerce (AREA)
- Complex Calculations (AREA)
Abstract
The invention discloses a singular optimal control simultaneous solving method based on partial moving finite element nodes. The method is characterized in that the switching points between the singular arcs and the nonsingular arcs are not easy to find accurately, the singular optimal control problem is dispersed into a nonlinear programming problem by using finite element orthogonal configuration, the positions of finite element nodes needing to be moved are determined by a switching function based on a control equation, and the positions of the finite element nodes needing to be moved are used as optimization variables and are optimized together with the control variables, so that the switching points between the singular arcs and the nonsingular arcs are found with high accuracy. For the singular optimal control problems existing in a large number in the production process, the method can obtain the high-precision optimal control curve, is high in calculation speed, and is suitable for solving the large-scale singular optimal control problems.
Description
Technical Field
The invention relates to the field of process system control, in particular to a simultaneous solving method based on partial moving finite element nodes.
Background
Singular optimal control problems are common in chemical, pharmaceutical and like manufacturing processes, such as catalyst mixing, metabolite optimization in fed-batch cultures, and the like. For the Singular Optimal Control problem frequently appearing in process system Control, a classical processing method is a regularization method, Jacobson proposes to add penalty factors in an objective function and changes the solving of the Singular Optimal Control problem into the solving of a series of nonsingular Optimal Control subproblems [ Jacobson D, Gershwin S, Lele M.computation of Optimal simple Controls [ J ]. IEEE Transactions on Automatic Controls, 1970,15(1):67-73 ]. However, this method is computationally cumbersome, computationally expensive, and results in numerical difficulties when the penalty factor approaches 0. Even if the sub-problems of the regularization method are solved by using the existing commercialized pseudo-spectrum method, a good result cannot be obtained. In addition, in order to improve the calculation precision, the pseudo-spectral Method adds a large number of grid nodes or increases the order of an interpolation polynomial in the Solving process, so the calculation amount for Solving the sub-problem is also large [ Darby CL, Hager WW, Rao AV. an hp-Adaptive pseudo-spectral Method for Solving Optimal Control schemes [ J ]. Optimal Control Applications and Methods,2011,32(4):476-502 ].
Disclosure of Invention
The invention aims to provide a singular optimal control simultaneous solving method based on partial moving finite element nodes aiming at the defects of the prior art.
The purpose of the invention is realized by the following technical scheme: a singular optimal control simultaneous solving method based on partial moving finite element nodes comprises the following steps:
the method comprises the following steps: equally dividing a control time domain into M finite elements, adopting P orthogonal configuration points on each finite element, dispersing an original singular optimal control problem into a nonlinear programming problem, and solving the nonlinear programming problem to obtain an approximate optimal control curve;
step two: detecting the discrete error of the non-configuration points between the adjacent configuration points on each finite element, if the discrete error is less than the specified tolerance epsilon1Then go on step three; otherwise, inserting finite element nodes on the non-configuration points, and re-solving the nonlinear programming problem until the discrete error on the non-configuration points meets the tolerance epsilon1Until now, the final finite element number obtained in the step is recorded as N, and the finite element node distribution is recorded as a1，a2，…，aN；
Step three: calculating a switching function on each finite element based on the obtained approximate optimal control curve, inserting a new finite element node if the value of the switching function is not accordant with the theoretical value, and optimizing the newly inserted finite element node and the control variable which are used as variables to be optimized together to obtain an improved optimal control curve;
step four: based on the improved optimal control curve, if the slope of the control curve on the jth finite element is larger than a preset threshold value, deleting the finite element, otherwise, finishing the calculation;
step five: if the finite element is deleted in the fourth step, fixing all the remaining finite element nodes, re-optimizing the control curve, and finishing the calculation.
Further, the third step is realized by the following sub-steps:
(3.1) if the switching function is not zero and the control variables are not at the upper and lower boundaries on the ith finite element, then in that finite element [ a ]i-1,ai]Inserting new finite element node b in the middle pointi；
(3.2) fixing the original node a0,a1,…,aNIf there are K finite element nodes (marked as K) newly inserted
) Position of the node
As optimization variables, with control variables, and requires
Satisfy the requirement of
Compared with the traditional method for processing the singular optimal control problem, the method has the advantages that the switching points between the singular arcs and the nonsingular arcs can be found more accurately, so that a better optimal control curve can be found, in addition, the calculation speed is high, and the method is suitable for solving the large-scale singular optimal control problem.
Drawings
FIG. 1 is a flow chart of a singular optimal control simultaneous solution method based on partially moved finite element nodes;
FIG. 2 is an optimal control curve obtained by solving a catalyst mixing case using a singular optimal control simultaneous solution method based on partially moving finite element nodes;
FIG. 3 is an optimal control curve obtained by solving a catalyst mixing case by a regularization and pseudo-spectrum method;
FIG. 4 is an optimal control curve obtained by solving a Lee-Ramirez bioreactor case using a singular optimal control simultaneous solution method based on partially moving finite element nodes;
FIG. 5 is an optimal control curve obtained by solving the Lee-Ramirez bioreactor case by using a regularization and pseudo-spectral method.
Detailed Description
The present invention is described in detail below with reference to the accompanying drawings.
As shown in fig. 1, the singular optimal control simultaneous solution method based on partial moving finite element nodes of the present invention includes the following steps:
the method comprises the following steps: equally dividing a control time domain into M finite elements, adopting P orthogonal configuration points on each finite element, dispersing an original singular optimal control problem into a nonlinear programming problem, and solving the nonlinear programming problem to obtain an approximate optimal control curve;
step two: detecting the discrete error of the non-configuration points between the adjacent configuration points on each finite element, if the discrete error is less than the specified tolerance epsilon1Then go on step three; otherwise, inserting finite element nodes on the non-configuration points, and re-solving the nonlinear programming problem until the discrete error on the non-configuration points meets the tolerance epsilon1Until now, the final finite element number obtained in the step is recorded as N, and the finite element node distribution is recorded as a1，a2，…，aN；
Step three: and calculating a switching function on each finite element based on the obtained approximate optimal control curve, inserting a new finite element node if the value of the switching function is not accordant with the theoretical value, and optimizing the newly inserted finite element node and the control variable which are taken as variables to be optimized together to obtain the improved optimal control curve.
The step is the core of the invention and is divided into the following substeps:
1) if the absolute value of the switching function on the ith finite element is larger than the tolerance epsilon2And the distance between the control variable and the upper boundary and the distance between the control variable and the lower boundary are both larger than the tolerance epsilon3Then in the finite element [ a ]i-1,ai]Inserting new finite element node b in the middle pointi；
2) Fixing the original node a0,a1,…,aNIf there are K finite element nodes (marked as K) newly inserted
) Position of the node
As optimization variables, with control variables, and requires
Satisfy the requirement of
Step four: based on the improved optimal control curve, if the slope of the control curve on the jth finite element is larger than the predetermined threshold epsilon4If not, deleting the finite element, otherwise, finishing the calculation;
step five: if the finite element is deleted in the fourth step, fixing all the remaining finite element nodes, re-optimizing the control curve, and finishing the calculation.
Examples
In the following two embodiments, the parameters involved in the invention are set as: m is 9, P is 5, epsilon1＝10-5，ε2＝10-6，ε3＝10-3，ε4＝100。
Example 1
The following reaction is assumed:
a and B are reversibly catalyzed by catalyst I, and B to C are irreversibly catalyzed by catalyst II. How we need to consider mixing the two catalysts to maximize the yield of product C can be described as the proposition:
max1-a(tf)-b(tf)
a(0)＝1,b(0)＝0,u(t)∈[0,1]
wherein tf represents the total length of the reactor, t represents the distance to the reactor inlet, a (t) and B (t) represent the mole fractions of reactants A and B,
and
represents the rate of change to t; k1, k2 and k3 represent reaction rate constants for reactions 1,2 and 3; u (t) represents the proportion of catalyst type I.
In this example, tf＝4,k1＝k3＝1,k 210, comprising the following steps:
the method comprises the following steps: dividing the control time domain into 9 finite elements equally, wherein the corresponding finite element nodes are as follows: 0.44444, 0.88889, 1.33333, 1.77778, 2.22222, 2.66667, 3.11111, 3.55556 and 4.00, dispersing the original singular optimal control problem into a nonlinear programming problem by adopting an orthogonal configuration method, and solving the nonlinear programming problem to obtain an approximate optimal control curve;
step two: detecting a dispersion error at a non-configuration point between adjacent configuration points on each finite element, requiring less than a specified tolerance epsilon1The number of finite elements finally obtained in the step is 15;
step three: calculating a switching function on each finite element based on the obtained approximate optimal control curve, inserting a new finite element node if the value of the switching function is not accordant with the theoretical value, and optimizing the newly inserted finite element node and the control variable which are taken as variables to be optimized together to obtain the improved optimal control curve:
1) the absolute value of the switching function in the 3 rd, 4 th, 5 th, 11 th, 12 th, 13 th and 14 th finite elements is larger than epsilon2And the distance between the control variable and the upper boundary and the lower boundary is larger than epsilon3Then in finite elements [0.074211,0.19121 ]]，[0.19121,0.320851]，[0.320851,0.413386]，[2.2222,2.66664]，[2.66664,3.111084]，[3.111084,3.55552]，[3.55552,3.99996]Inserting a new finite element node into the middle point;
2) fixing the original finite element nodes, optimizing the newly inserted 7 finite element nodes 0.1363, 0.25603, 0.3672, 2.44442, 2.88886, 3.3333 and 3.72523 serving as optimization variables together with control variables, and requiring the newly inserted finite element nodes to be always in 7 intervals listed in 1);
step four: based on the improved optimal control curve, the control curve slopes on all finite elements are found to be smaller than a preset threshold value 100, the calculation is finished, the obtained final control curve is shown in fig. 2, and the switching point between the singular arc and the non-singular arc in the result is almost completely coincided with the theoretical value.
In order to further illustrate the beneficial effects of the invention, a method based on the combination of regularization and pseudo-spectrum method is also adopted to solve the problem, and the result is shown in fig. 3 and is far from the theoretical result.
Example 2
This example is an optimal control problem for a Lee-Ramirez bioreactor, with the goal of maximizing the amount of product by controlling the glucose and inducer feed rates, and the specific examples are as follows:
min-x4(tf)x1(tf)
s.t.
t3(t)＝0.22x7(t)/t2(t)+x6(t)
g1(t)＝x3(t)t3(t)/t1(t)；
g2(t)＝0.233x3(t)[0.0005+x5(t)]/[0.022+x5(t)]/t1(t)；
g3(t)＝0.09x5(t)/[0.034+x5(t)]；0≤u1(t),u2(t)≤1
wherein t represents the reaction time, tf represents the total reaction time, and x1(t) represents the volume of the reactor, x2(t) represents cell density, x3(t) represents the nutrient concentration, x4(t) represents the concentration of a foreign protein, x5(t) denotes inducer concentration, x6(t) denotes the inducer shock factor, x7(t) means an inducer recovery factor; u. of1(t) represents the glucose feed rate, u2(t) represents inducer feed flow; c. C1Denotes the feed concentration of glucose, c2Represents the coefficient of increase in yield, c3Represents the feed concentration of the inducer; g1(t) represents the specific growth rate, g2(t) shows the productivity of the foreign protein, g3(t) denotes shock rate; t is t1(t)，t2(t) and t3(t) is an intermediate variable.
The method comprises the following steps:
the method comprises the following steps: dividing a control time domain into 9 finite elements, wherein the corresponding finite element nodes are as follows: 0,1.11111, 2.22222, 3.33333, 4.44444, 5.55556, 6.66667, 7.77778, 8.88889 and 10, dispersing the original singular optimal control problem into a nonlinear programming problem by adopting an orthogonal configuration method, and solving the nonlinear programming problem to obtain an approximate optimal control curve;
step two: detecting a dispersion error at a non-configuration point between adjacent configuration points on each finite element, requiring less than a specified tolerance epsilon1The number of finite elements finally obtained in the step is 35;
step three: calculating a switching function on each finite element based on the obtained approximate optimal control curve, inserting a new finite element node if the value of the switching function is not accordant with the theoretical value, and optimizing the newly inserted finite element node and the control variable which are taken as variables to be optimized together to obtain the improved optimal control curve:
1) the absolute value of the switching function in 1,2,3,4,15,16,17,18,19,26,27 finite elements is greater than epsilon2And the distance between the control variable and the upper boundary and the lower boundary is larger than epsilon3Then in finite element [0,1.11111 ]]，[1.11111，2.22222]，[2.22222,3.33333]，[3.33333,4.44444]，[4.72385,5.55558]，[5.55558,5.58731]，[5.58731,6.66671]，[6.66671,7.77781]，[7.77781,7.809533]，[8.26496,8.353871]，[8.35387,8.7103]Inserting a new finite element node into the middle point;
2) fixing the original finite element nodes, optimizing 11 newly inserted finite element nodes 1.10967, 1.11252, 2.22367, 3.33530, 4.82270, 5.55662, 6.66563, 7.26963, 7.77908, 8.3526 and 8.6585 serving as optimization variables and control variables, and requiring the newly inserted finite element nodes to be always in 11 intervals listed in 1);
step four: based on the improved optimal control curve, if the slope of the control curve on the 2 nd, 3 rd, 7 th, 20 th, 24 th, 26 th and 36 th finite elements is larger than a preset threshold value 100, deleting the finite elements;
step five: if the finite element is deleted in the fourth step, all the remaining finite element nodes are fixed to 39, the control curve is re-optimized, the calculation is ended, the final control curve is shown in fig. 4, and the final target value is-6.1516968538.
To further illustrate the beneficial effects of the invention, a method based on regularization combined with pseudo-spectrometry was also used to solve the problem, with the result shown in FIG. 5, where the target value is-6.1515979885. The objective function value is comparatively poor.
Claims (2)
1. A singular optimal control simultaneous solving method based on partial moving finite element nodes is characterized by comprising the following steps:
the method comprises the following steps: equally dividing a control time domain into M finite elements, adopting P orthogonal configuration points on each finite element, dispersing an original singular optimal control problem into a nonlinear programming problem, and solving the nonlinear programming problem to obtain an approximate optimal control curve;
step two: detecting the discrete error of the non-configuration points between the adjacent configuration points on each finite element, if the discrete error is less than the specified tolerance epsilon1Then go on step three; otherwise, inserting finite element nodes on the non-configuration points, and re-solving the nonlinear programming problem until the discrete error on the non-configuration points meets the tolerance epsilon1Until now, the final finite element number obtained in the step is recorded as N, and the finite element node distribution is recorded as a1，a2，…，aN；
Step three: calculating a switching function on each finite element based on the obtained approximate optimal control curve, inserting a new finite element node if the value of the switching function is not accordant with the theoretical value, and optimizing the newly inserted finite element node and the control variable which are used as variables to be optimized together to obtain an improved optimal control curve;
step four: based on the improved optimal control curve, if the slope of the control curve on the jth finite element is larger than a preset threshold value, deleting the finite element, otherwise, finishing the calculation;
step five: and if the finite element is deleted in the fourth step, fixing all the remaining finite element nodes, re-optimizing the control curve, and finishing the calculation.
2. The singular optimal control simultaneous solution method based on partially moved finite element nodes according to claim 1, wherein the step three is realized by the following sub-steps:
(3.1) if the absolute value of the switching function in the ith finite element is larger than the tolerance ε2And the distance between the control variable and the upper boundary and the distance between the control variable and the lower boundary are both larger than the tolerance epsilon3Then in the finite element [ a ]i-1,ai]Inserting new finite element node b in the middle pointi；
(3.2) fixing the original node a0,a1,…,aNIf there are K finite element nodes (marked as K) newly inserted
Node position
As optimization variables, with control variables, and requires
Satisfy the requirement of
Priority Applications (1)
| Application Number | Priority Date | Filing Date | Title |
|---|---|---|---|
| CN201910120094.9A CN109814391B (en) | 2019-02-18 | 2019-02-18 | Singular optimal control simultaneous solving method based on partial moving finite element nodes |
Applications Claiming Priority (1)
| Application Number | Priority Date | Filing Date | Title |
|---|---|---|---|
| CN201910120094.9A CN109814391B (en) | 2019-02-18 | 2019-02-18 | Singular optimal control simultaneous solving method based on partial moving finite element nodes |
Publications (2)
| Publication Number | Publication Date |
|---|---|
| CN109814391A CN109814391A (en) | 2019-05-28 |
| CN109814391B true CN109814391B (en) | 2021-09-17 |
Family
ID=66606834
Family Applications (1)
| Application Number | Title | Priority Date | Filing Date |
|---|---|---|---|
| CN201910120094.9A Active CN109814391B (en) | 2019-02-18 | 2019-02-18 | Singular optimal control simultaneous solving method based on partial moving finite element nodes |
Country Status (1)
| Country | Link |
|---|---|
| CN (1) | CN109814391B (en) |
Citations (4)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| CN104250034A (en) * | 2014-03-24 | 2014-12-31 | 杭州电子科技大学 | Operation optimization method of full flow roll type reverse osmosis seawater desalination system |
| CN104573236A (en) * | 2015-01-08 | 2015-04-29 | 浙江大学 | Continuous system simulation verification-based partial finite element method |
| CN104696944A (en) * | 2015-01-26 | 2015-06-10 | 浙江大学 | Dynamic optimization and parameter estimation integrated method based on load prediction |
| CN105335797A (en) * | 2015-11-03 | 2016-02-17 | 浙江大学 | Automatic parking locus optimization method based on full-simultaneous dynamic optimization framework |
-
2019
- 2019-02-18 CN CN201910120094.9A patent/CN109814391B/en active Active
Patent Citations (4)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| CN104250034A (en) * | 2014-03-24 | 2014-12-31 | 杭州电子科技大学 | Operation optimization method of full flow roll type reverse osmosis seawater desalination system |
| CN104573236A (en) * | 2015-01-08 | 2015-04-29 | 浙江大学 | Continuous system simulation verification-based partial finite element method |
| CN104696944A (en) * | 2015-01-26 | 2015-06-10 | 浙江大学 | Dynamic optimization and parameter estimation integrated method based on load prediction |
| CN105335797A (en) * | 2015-11-03 | 2016-02-17 | 浙江大学 | Automatic parking locus optimization method based on full-simultaneous dynamic optimization framework |
Non-Patent Citations (1)
| Title |
|---|
| Trajectory optimization for lunar soft landing with a;Lin Ma等;《Advances in Engineering Software》;20160801;第265-276页 * |
Also Published As
| Publication number | Publication date |
|---|---|
| CN109814391A (en) | 2019-05-28 |
Similar Documents
| Publication | Publication Date | Title |
|---|---|---|
| Lee et al. | Optimal fed‐batch control of induced foreign protein production by recombinant bacteria | |
| CN109696827A (en) | The pid parameter setting method of inertia weight cosine adjustment particle swarm optimization algorithm | |
| CN110426956B (en) | Intermittent process optimal compensation control strategy based on process migration model | |
| Peroni et al. | Optimal control of a fed-batch bioreactor using simulation-based approximate dynamic programming | |
| CN103088448B (en) | A kind of carbon fibre precursor jet stretch technique controlled based on data-driven cooperative intelligent | |
| CN108919639A (en) | A kind of PID controller parameter best proportion method for establishing model | |
| CN105353607A (en) | Data-difference-driving batch process self-learning dynamic optimization method | |
| CN109814391B (en) | Singular optimal control simultaneous solving method based on partial moving finite element nodes | |
| CN111123708B (en) | Coking furnace hearth pressure control method based on distributed dynamic matrix control optimization | |
| CN112398170A (en) | Parameter identification method for photovoltaic inverter system model | |
| CN103412486A (en) | Intelligent control method for polyvinyl chloride steam stripping process | |
| CN111610751B (en) | Interpolation error multi-subdivision iterative calculation method for cross point set NURBS interpolation curve | |
| CN108388218A (en) | The adaptive batch process optimization method of amendment based on latent variable process migration model | |
| CN111506037A (en) | Dynamic matrix optimization distributed control method for industrial heating furnace system | |
| CN109635330B (en) | Method for accurately and rapidly solving complex optimization control problem based on direct method | |
| CN111427261A (en) | PID parameter setting method based on cat swarm algorithm | |
| CN115857432A (en) | Cutter track smoothing method and system based on curvature maximum optimization | |
| Würth et al. | On the numerical solution of discounted economic NMPC on infinite horizons | |
| CN110109430A (en) | A kind of intermittent beer fermenting device Optimal Control System | |
| CN116149185A (en) | Double undisturbed switching control method of sewage treatment system | |
| Lu et al. | Stability and fuel economy of nonlinear vehicle platoons: A distributed economic MPC approach | |
| Blanch et al. | Optimal conditions for gramicidin S production in continuous culture | |
| CN115079569A (en) | AGV nonlinear active disturbance rejection control method based on longicorn stigma search algorithm | |
| Ma | A linear optimal feedback control for producing 1, 3-propanediol via microbial fermentation | |
| CN104298213A (en) | Index time varying gain type iterative learning control algorithm based on reference batch |
Legal Events
| Date | Code | Title | Description |
|---|---|---|---|
| PB01 | Publication | ||
| PB01 | Publication | ||
| SE01 | Entry into force of request for substantive examination | ||
| SE01 | Entry into force of request for substantive examination | ||
| GR01 | Patent grant | ||
| GR01 | Patent grant |