- Email: [email protected]
Flatness-Based Model Predictive Control with Linear Programming for a Single Mast Stacker Crane
Proceedings of the 9th Vienna International Conference on Proceedings of the 9th Vienna International Conference on Mathematical Proceedings ofModelling the 9th Vienna International Conference on Proceedings ofModelling the 9th Vienna International Conference on Mathematical Vienna, Austria, February 21-23, 2018 Available online at www.sciencedirect.com Mathematical Modelling Mathematical Modelling Vienna, Austria, February 21-23, 2018 Vienna, Austria, February 21-23, 2018 Vienna, Austria, February 21-23, 2018
ScienceDirect
IFAC PapersOnLine 51-2 (2018) 31–36
Flatness-Based Model Predictive Control Flatness-Based Model Predictive Control Flatness-Based Model Predictive Control Flatness-Based Model Predictive Control with Linear Programming for a Single Mast with Linear Programming for a Single Mast with Linear Programming for a Single Mast with Linear Programming for a Single Mast Stacker Crane Stacker Crane Stacker Crane Stacker Crane ∗ ∗
Galkina, A. ∗ Schlacher, K. ∗ Galkina, A. ∗∗ Schlacher, K. ∗∗ Galkina, Galkina, A. A. Schlacher, Schlacher, K. K. ∗ ∗ Institute of Automatic Control and Control Systems Technology, ∗ Institute of Automatic Control and Control Systems Technology, of Control Control Systems Technology, Johannes University Linz, and Altenberger 4040 Linz, ∗ Institute InstituteKepler of Automatic Automatic Control and Control Strasse Systems69, Technology, Johannes Kepler University Linz, Altenberger Strasse 69, 4040 Linz, Johannes Kepler University Linz, Altenberger Strasse 69, Austria,(e-mail: [email protected]). Johannes Kepler University Linz, Altenberger Strasse 69, 4040 4040 Linz, Linz, Austria,(e-mail: [email protected]). Austria,(e-mail: [email protected]). Austria,(e-mail: [email protected]). Abstract: This paper deals with a trajectory tracking model predictive control for a single Abstract: paper deals with a trajectory model predictive control for aa single Abstract: This paper deals with a tracking model predictive control for mast stackerThis crane, which is used automatic tracking storage or retrieval of payloads automated Abstract: This paper deals with for a trajectory trajectory tracking model predictive controlin for a single single mast stacker crane, which is used for automatic storage or retrieval of payloads in automated mast stacker crane, which is used for automatic storage or retrieval of payloads in automated warehouses. The mathematical model of the plant is a distributed parameter one, but it admits mast stackerThe crane, which is used for of automatic storage or retrieval of payloads in automated warehouses. mathematical model the plant is a distributed parameter one, but it admits warehouses. mathematical the is one, it an excellent The approximation by model a flat of lumped parameter system. parameter The approximate nonlinear warehouses. The mathematical model of the plant plant is aa distributed distributed parameter one, but but nonlinear it admits admits an excellent approximation by a flat lumped parameter system. The approximate an excellent approximation by a flat lumped parameter system. The approximate nonlinear system is simplified further and the resulting linear time-varying system is implemented to an excellent approximation by a flat lumped parameter system. The approximate nonlinear system is simplified further and the resulting linear time-varying system is implemented to system is simplified further and the resulting linear time-varying system is implemented to design a model predictive controller. To provide an efficient formulation of an optimization task system ismodel simplified further and theToresulting linear time-varying system implementedtask to design predictive controller. provide an efficient formulation of anisoptimization design a model controller. To provide formulation of task for thea the LTV-model parametrized by a special flat output. Then, the model design acontroller, model predictive predictive controller.is To provide an an efficient efficient formulation of an an optimization optimization task for the controller, the LTV-model is parametrized by a special flat output. Then, the model for the LTV-model is by special output. Then, the predictive control isthe formulated in the form of a linear For solving derived linear for the controller, controller, the LTV-model is parametrized parametrized by aaprogram. special flat flat output. the Then, the model model predictive control is formulated in the form of aa linear program. For solving the derived linear predictive control is formulated in the form of linear program. For solving the derived linear program, an open-source optimal solver lp solve is chosen. To reduce a mechanical stress at predictive control is formulated in solver the form of a linear program. For solving the derived linear program, an open-source optimal lp solve is chosen. To reduce a mechanical stress at program, an open-source optimal solver lp solve is chosen. To reduce a mechanical stress clamping of a mast and an overall accessibility time, a special optimal trajectory calculated program, an open-source optimal solver lp solvetime, is chosen. To reduce a trajectory mechanicalcalculated stress at at clamping of a mast and an overall accessibility aa special optimal clamping of a mast and an overall accessibility time, special optimal trajectory calculated offline is chosen for the tracking by the designed controller. Finally, the simulation results of clamping of a mast and tracking an overall accessibility time, a special optimal trajectory calculated offline is chosen for the by the designed controller. Finally, the simulation results of offline is chosen for the tracking by the designed controller. Finally, the simulation results the nonlinear approximate model, subjected to model uncertainties and external disturbances, offline is chosen for the tracking by the designed controller. Finally,and theexternal simulation results of of the nonlinear approximate model, subjected to model uncertainties disturbances, the nonlinear approximate model, subjected model uncertainties and external disturbances, with the designed model predictive controllerto are presented. the nonlinear approximate model, subjected to model uncertainties and external disturbances, with the designed model predictive controller presented. with the designed model controller are are presented. with theIFAC designed model predictive predictive are presented. © 2018, (International Federation controller of Automatic Control) Hosting by Elsevier Ltd. All rights reserved. Keywords: Model predictive control, Flatness, Linear programming, Flexible construction Keywords: Model predictive control, Flatness, Linear programming, Flexible construction Keywords: Keywords: Model Model predictive predictive control, control, Flatness, Flatness, Linear Linear programming, programming, Flexible Flexible construction construction 1. INTRODUCTION of modeling, see e.g. Bachmayer et al. (2011), considers the 1. INTRODUCTION INTRODUCTION of modeling, modeling, see e.g.during Bachmayer et al. al. (2011), considersThe the 1. of e.g. Bachmayer et considers the payload beingsee fixed positioning of the carriage. 1. INTRODUCTION of modeling, see e.g.during Bachmayer et al. (2011), (2011), considersThe the payload being fixed positioning of the carriage. being during positioning of the carriage. The second type of fixed the modeling approach assumes the payload A flexible mechanical system, depicted in Fig.1, represents payload payload being fixed during positioning of the carriage. The secondmoved type of of the the its modeling approach assumes the payload payload flexiblemast mechanical system, depicted in type modeling assumes the being positionapproach is treated as an uncertainty. A flexible mechanical depictedalso in Fig.1, Fig.1, represents aA single stacker system, crane (SMC), called represents a flexible second second type ofand the its modeling approach assumes the payload being moved and position is treated as an uncertainty. A flexiblemast mechanical system, depicted in Fig.1, represents a single stacker crane (SMC), also called a flexible being moved and its position is treated as an uncertainty. In Hajdu andand Gaspar (2014) this approach with arack single mast stacker crane also called flexible feeder. is usually used for automatic being moved its position is treated as antogether uncertainty. Hajdu and modeling Gaspar (2014) this approach together with a single mastSMC stacker crane (SMC), (SMC), also called aa storage flexible In rack feeder. SMC is usually usually used for for automatic storage Hajdu and Gaspar (2014) this approach together with aIn multi-body technique is presented. The third rack feeder. SMC is used automatic storage or retrieval of payloads in automated warehouses. The In Hajdu and Gaspar (2014) this approach together with a multi-body modeling technique is presented. The third rack feeder. SMC is usually used for automatic storage or retrieval retrieval of payloads payloads in initial automated warehouses. The atype multi-body modeling technique is presented. The third third of the modeling approach exploiting a fast underlying or of in automated warehouses. The payload is moved from the position to the desired a multi-body modeling technique is presented. The type of of the the modeling approach exploiting fast underlying or retrieval of payloads in initial automated warehouses. The type payload is moved moved from the the position toreduced the desired desired modeling exploiting aaa fast underlying control on approach the current converter both for the payload is from position to the end position. The weight of initial the mast is often due velocity type of the modeling approach exploiting fast underlying velocity control on the current converter both for and the payload is moved from the initial position toreduced the desired end position. The weight of the mast is often due velocity control on the current converter both for the payload and carriage cancurrent be found in Aschemann end position. The weight of the mast is often reduced due to economic reasons. Therefore, vibrations are generated velocity control on the converter both for and the payload and carriage can be found in Aschemann end position. The weight of the mast is often reduced due to economic economic reasons. Therefore, Therefore, vibrations areHence, generated payload and carriage can be found in Aschemann and Schindele (2010). Another modeling approach is presented to reasons. vibrations are generated during transportation and positioning process. the payload and carriage can be found in Aschemann and Schindele (2010). Another modeling approachassumes is presented presented to economic reasons. Therefore, vibrations areHence, generated during transportation and positioning positioning process. (2010). modeling approach is in Staudecker et Another al. (2008). This approach the during transportation and process. Hence, the Schindele structural vibrations must be taking into account in the Schindele (2010). Another modeling approachassumes is presented in Staudecker et al. (2008). This approach the during transportation and positioning process. Hence, the structural vibrations must be taking into account in the in Staudecker et al. (2008). This approach assumes the controlled motion of the payload during the carriage transstructural vibrations must be taking into account in the model being derived for control of the positioning process. in Staudecker et al. (2008). This approach assumes the controlled motion of the the payload during the carriage carriage transstructural vibrations must be taking into accountprocess. in the controlled model being derived for control of the positioning motion of payload during the transportation. The first and second approaches do not allow model control of the positioning process. SMC isbeing an derived essentialfor element of automated warehouses controlled motion of the payload during the carriage transportation. The first and second approaches do not allow model being derived for control of the positioning process. SMCcontributes is an an essential element of automated automated warehouses The first second not minimizing the overall accessibility time of thedo storage and SMC is of warehouses and a lotelement to its overall accessibility time. portation. portation. The first and and second approaches approaches not allow allow minimizing the overall overall accessibility time of of the thedostorage storage and SMC is an essential essential element of automated warehouses and contributes a lot to its overall accessibility time. minimizing the accessibility time and retrieval processes. and contributes aa lot to overall accessibility time. Therefore, the control to position the system as minimizing the overall accessibility time of the storage and retrieval processes. and contributes lotpurpose to its its is overall accessibility time. Therefore, the control purpose is to position the system as processes. Therefore, the purpose is position system as retrieval fast as possible. To achieve this purpose thethe feedforward retrieval processes. This paper deals with MPC for stabilization of the offline Therefore, the control control purpose is to to position system as fast as as possible. possible. To achieve achieve this purpose thethe feedforward This paper deals with with MPC for for stabilization ofthe themodel, offline fast To this purpose the feedforward control design based on modern control concepts such as This paper deals MPC stabilization the offline generated optimal trajectory and is based onof fast as possible. To achieve this purpose the feedforward control design design optimal based on ontrajectory modern control control concepts such as as This paperoptimal deals with MPC for stabilization ofthe themodel, offline generated trajectory and is based on control based modern concepts such flatness-based planning can be applied. optimal trajectory and on the model, which is presented in Staudecker etis al.based (2008). The mathecontrol design optimal based ontrajectory modern control concepts such as generated flatness-based planning can be applied. generated optimal trajectory and is based on the model, which is presented in Staudecker et al. (2008). The matheflatness-based optimal trajectory be applied. applied. But the considered system iscan subjected to maticalis model presented in Staudecker et The of the is a distributed parameter one, flatness-based optimalmechanical trajectory planning planning can be But the considered considered mechanical system is to which which is model presented in plant Staudecker et al. al. (2008). (2008). The mathemathematical of excellent the plant is aa distributed distributed parameter one, But the mechanical system is subjected subjected to parameter uncertainties and external disturbances such as matical model of the plant is parameter one, but it admits an approximation by a flat lumped But the considered mechanical system is subjected to parameter uncertainties and external disturbances such as as matical modelan of excellent the plantapproximation is a distributedbyparameter one, but it admits a flat lumped parameter uncertainties and external disturbances such friction force. Therefore, feedback extension is desired, to but it admits an excellent approximation by aa flat lumped parameter system. To reduce the mechanical stress due parameter uncertainties and external disturbances such as friction force. Therefore, feedback extension is desired, to but it admits an excellent approximation by flat lumped parameter system. To reduce the mechanical stress due friction force. Therefore, feedback extension is to stabilize the planned offline trajectory. Stabilization of the system. To reduce the mechanical stress due to a bending moment acting on the mast clamping and friction force. Therefore, feedback extension is desired, desired, to parameter stabilize the planned offline trajectory. Stabilization of the parameter system. To reduce the mechanical stress due to a bending moment acting on the mast clamping and stabilize the offline Stabilization of desired trajectory is the maintrajectory. objective of the contribution. aareduce bending acting on clamping and the moment transporting the mast optimal trajectory stabilize the planned planned offline Stabilization of the the to desired trajectory isthe the maintrajectory. objective of the contribution. contribution. to bending moment actingtime, on the the mast clamping and to reduce the transporting time, the optimal trajectory desired trajectory is the main objective of the For this purpose, design of a flatness-based model to reduce the transporting time, optimal described et al. (2017) is chosen to trajectory be stabidesired isthe the design main objective of the contribution. For this thistrajectory purpose, of linear flatness-based model to reduce by the Rams transporting time, the the optimal trajectory described by Rams etminimizes al. (2017) (2017) is bending chosen tomoment be stabistabiFor purpose, design of aa flatness-based model predictive control the (MPC) with programming is described by Rams et al. is chosen to be lized. This trajectory the at For this purpose, the design of a flatness-based model predictive control control (MPC) (MPC) with with linear linear programming programming is is described by Rams etminimizes al. (2017)the is bending chosen tomoment be stabilized. This trajectory at predictive proposed. lized. This trajectory minimizes the bending moment at clamping of the mast and, at the same time, moves predictive control (MPC) with linear programming is the proposed. lized. This trajectory minimizes the bending moment at the clamping of the mast and, at the same time, moves proposed. clamping of mast and, the time, moves payload from starting position to the poproposed. There may be found several approaches both for mathe- the the clamping of the thethe mast and, at at the same same time,end moves the payload from the starting position to the end poThere may may be found found several approaches both both for mathepayload the to position as fastfrom as possible. It is position worth mentioning, There be approaches mathematical modeling and several for feedforward/feedback control of the the payload the starting starting to the the end end that position as fast fastfrom as this possible. It is position worth mentioning, that There may be found several approaches both for for mathematical modeling and for feedforward/feedback control of sition as as possible. It is worth mentioning, that the presented in paper MPC concept can be applied matical modeling and for feedforward/feedback control of SMC in the literature. The modeling approaches differ in sition as fast as possible. It is worth mentioning, that the presented in this paper MPC concept can be applied matical and for feedforward/feedback control of the SMC in modeling the payload literature. The modeling approaches differ in presented in paper MPC concept can be applied with any trajectory could be suitable SMC in the literature. approaches differ the way andThe the modeling mast are treated, as well as in the in this thisand paper MPC concept for canthe be controlapplied withpresented any trajectory trajectory and could be suitable suitable for the controlSMC in the the payload literature. The modeling approaches differ in with the way the and the mast are treated, as well as in any and could be for the controlthe way the payload and the mast are treated, as well as in with any trajectory and could be suitable for the controlphysical variables being used as control inputs. One type the way the payload and the mast are treated, as well as in physical variables variables being used used as control control inputs. One One type physical physical variables being being used as as control inputs. inputs. One type type
Copyright © 2018, 2018 IFAC 1 Hosting by Elsevier Ltd. All rights reserved. 2405-8963 © IFAC (International Federation of Automatic Control) Copyright © 2018 IFAC 1 Copyright ©under 2018 responsibility IFAC 1 Control. Peer review of International Federation of Automatic Copyright © 2018 IFAC 1 10.1016/j.ifacol.2018.03.006
Proceedings of the 9th MATHMOD 32 Vienna, Austria, February 21-23, 2018
xk
A. Galkina et al. / IFAC PapersOnLine 51-2 (2018) 31–36
point mass mh moves along the beam centerline. The forces Fx and Fy acting on the carriage and the lifting unit, respectively, serve as control inputs.
mk
w(Y, t) xh
yh
Concerning a mathematical notation, the following shortcuts are used, namely ∂t = ∂/∂t and ∂Y = ∂/∂Y as well as ∂tk = ∂ k /∂tk and ∂Yk = ∂ k /∂Y k , k = 2, ..., n. When dealing with functions u = u(t), it is written u˙ = ∂t u, u ¨ = ∂t2 u and u(i) for the time derivatives of the higher order.
mh
Fy
Y X
Fx
2.1 System Modeling
mc
To shorten some expressions, we interchange sometimes the function, e.g. w (Y, t), and the coordinate, e.g. w, notation. Taking the mentioned above assumptions into account, the kinetic energy of the system is ˆL mk 2 ρA mc 2 mh 2 2 2 Ek = x˙ + x˙ h + y˙ h + x˙ + (∂t w) dY 2 c 2 2 k 2
xc Fig. 1. Single mast stacker crane oriented model presented, for example, in Aschemann and Schindele (2010) as well. The paper consists of six main parts. At first, to keep a consistency of the proposed approach, modeling of the infinite-dimensional system and the system approximation by a finite-dimensional model is recapitulated in Section 2. In general, the nonlinear approximate model can be used together with a nonlinear MPC for stabilization of a desired trajectory. But the optimization problem will be a nonlinear one, that increases the calculation time considerably and it is hardly applicable to a real-time system. For this reason, in Section 3, the flatness analysis, which is proposed in Staudecker et al. (2008), is implemented to simplify the approximate model. As a result, a linear timevarying (LTV) model is derived. We consider the problem of stabilization of the optimal trajectory being generated offline and solve it by means of model predictive control. For MPC design the simplified LTV-model must be discretized with respect to time. For this aim, discretization of a flat space is proposed in Section 4. Also, to improve an efficiency of MPC realized by the free open-source optimal solver lp solve, in Section 4 the discretized flat system states are parameterized by a flat output again. Afterward, in Section 5, the stabilizing MPC in the form of a linear program is presented. The flatness-based MPC design with linear programming exploiting the simplified system model of SMC represents the main contribution of this paper. Finally, the simulation results with the nonlinear finite-dimensional system, which is subjected to parameter uncertainties and external disturbances, are presented in Section 6.
0
and the potential energy is
1 Ep = mh gyh + EI 2
ˆL 0
∂Y2 w
2
dY.
The contribution of the external forces can be expressed as Eext = xc Fx + yh Fy . Hence, the Lagrangian of the system can be written as L = Ekin − Epot + Eext . To derive the equations of motion for the SMC the Hamilton’s principle is applied and the following equations of motion are derived in the form of one partial differential equation ρA∂t2 w = −EI∂Y4 w
and the set of ordinary differential equations mc x ¨c = Fx − EI∂Y3 w (0, t) mh y¨h = Fy − ¨h ∂Y w (y h , t) mh g − mh x mh x ¨h = EI ∂Y3 w yh− , t − ∂Y3 w yh+ , t ¨k = EI∂Y3 w (L, t) . mk x Also, the boundary conditions must be satisfied w (0, t) = xc , w (L, t) = xk , ∂Y w (0, t) = 0 , ∂Y2 w (L, t) = 0 , w (yh , t) = xh , ∂Y2 w yh− , t = ∂Y2 w yh+ , t .
(1)
(2)
(3)
For modeling details see Staudecker et al. (2008). The values yh+ and yh− denote the most upper and the most lower position yh of the payload. Next, the mixed-dimensional system equations (1), (2) and (3) are discretized in accordance with the Rayleigh-Ritz method and the finitedimensional approximation is derived.
2. MODELING In this section, we recapitulate briefly the derivation of the mathematical model, which is proposed in Staudecker et al. (2008). The carriage, the lifting unit and the tip with the corresponding masses mc , mh and mk are modeled as rigid bodies with the mechanical degrees of freedom xc , xh , yh and xk in the fixed coordinate frame XY (see Fig. 1). The displacement of the mast centerline is denoted by w (Y, t) and it is assumed that the mast can be considered as a beam, which satisfies the Euler-Bernoulli hypothesis, having length L, uniform mass density ρ, uniform cross section surface A and uniform flexural rigidity EI. The
2.2 System Approximation Applying the Rayleigh-Ritz method, a finite-dimensional approximation of the system model, represented by (1), (2) and (3), can be derived. According to this method the mast deflection w (Y, t) is approximated by the first-order Ritz ansatz function w∗ (Y, t) = xc (t) + Φ1 (Y ) q¯1 (t) (4) 2
Proceedings of the 9th MATHMOD Vienna, Austria, February 21-23, 2018
A. Galkina et al. / IFAC PapersOnLine 51-2 (2018) 31–36
only the function of yh,d . Substituting the simplification (9) into (7) yields the LTI-system 0 1 0 ˙y ¯= ¯+ y v¯ (10) 0 0 1
with the new generalized coordinate q¯1 , and the appropriately chosen spatial basis function 4
with the input v¯ = (Fy − gmh )/mh . Hence, the resulting simplified system can be written in the following form ¯˙ = Ax (yh,d ) x ¯ + bx (yh,d ) u x ¯ (11) ¯ + by v¯ ¯˙ = Ay y y with the LTV- and LTI-subsystems. The numerical verification (see Fig.(2)) shows, that applying the same input to the approximate (5) and simplified approximate model (11), the only small difference in the state dynamics occurs. Hence, the considered simplified model (11) can be used for the design of MPC in the form of a linear program. XY - payload position
0.6
Y (m)
Then, based on the flatness analysis, which is also described in Staudecker et al. (2008), the simplification of the approximate system (5) is carried out.
0.4 0.2 0
3. SYSTEM SIMPLIFICATION
0
0.2
0.4
0.6
0.8
Mb - bending moment
1
Exact AM Simplified AM
Mb (Nm)
3
Φ1 (Y ) = 6 (Y /L) − 4 (Y /L) + (Y /L) . Then, using the Ritz ansatz function (4) and applying the Hamilton’s principle, the approximate equations of motion can be written in the form of the explicit nonlinear ordinary differential equation of the second order ¨ + C (q, q) ˙ = Gu. M (q) q (5) The derived approximate system has the generalized co T ordinate q = xc q¯1 yh and the control input u = T [ Fx Fy ] . The bending moment acting at the clamping of the mast is proportional to q¯1 and can be expressed as Mb = 12EI q¯1 /L2 . The details and verification of the considered approximation can be found in Staudecker et al. (2008).
0.5 0 -0.5
1
-1
Exact AM Simplified AM
0
X (m) y˙ h - payload velocity
0.4 y˙ h (m/s)
In spite of the fact, that the state space representation of (5) is not state feedback equivalent to a linear timeinvariant (LTI) system, nonetheless it is flat. In Staudecker et al. (2008) it is shown, how the system decomposition can be used for the flatness analysis of the approximate system (5). It is also shown that yh is a flat output corresponding to the vertical dynamics and the desired trajectory yh,d (t) with the corresponding time derivatives are known (are computed offline). The system (5) is decomposed into a linear time-varying and nonlinear subsystems, namely ¯˙ = A (yh,d , y˙ h,d , y¨h,d ) x ¯ + b (yh,d ) u x ¯ (6) 1 1 T ˙ ¯ = xc q¯ x˙ c q¯ with the state vector x , the input u ¯ = Fx and ¯˙ = f (¯ y y, u ˆ) (7)
2
0.2 0 -0.2 -0.4
Exact AM Simplified AM
0
0.5
1
1.5 t (s)
2
˙ b (Nm/s) M
2
33
2.5
0.5 1 1.5 2 t (s) ˙ b - bending moment rate of change M
2.5
0 -2 -4
Exact AM Simplified AM
0
0.5
1
1.5
2
2.5
t (s)
Fig. 2. Approximate model (AM) (5) and simplified approximate model (11) 4. FLATNESS-BASED PARAMETRIZATION For the formulation of an optimization problem a discretetime system representation of (11) must be found. A dynamical system can be discretized by means of an exact or non-exact discretization method. To discretize the LTVsystem and to derive the efficient formulation of an optimization problem in the form of a linear program, the approach of a flatness-based parameterization of the system states, presented in Galkina et al. (2017) is used. We transform the simplified system into Brunovsky canonical form. In other words, the transformation of LTV-system into LTI-system is used. Then the integrator chain is discretized with respect to time by means of the exact time discretization and the discretized system is transformed into Brunovsky canonical form one more time.
¯ = [ yh y˙ h ] and the input u with the state y ˆ = T 1 1 xc q¯ x˙ c q¯˙ Fx Fy . The system (6) can be treated as a linear time-varying one and the corresponding methods and approaches for LTV-system can be implemented. Since for a linear system the properties of flatness and controllability are equivalent, we may obtain the second flat output corresponding to the dynamics in X-direction by transforming (6) into a time-varying controllability canonical form. Now, for the considered approximate model of SMC further simplifications can be made, to simplify the formulation of the linear optimization program for the stabilizing MPC. It is worth mentioning that the LTVtransformation 4 Tj1 (t) x ¯j (8) h1 = T
4.1 Motivation For the exact time discretization, a solution of a dynamic system is needed. In the case of LTI-system, a general solution is always available. Conversely, LTV- and nonlinear systems do not have a general solution. Both LTI- and LTV-systems need for a transition matrix to get a solution. In general, the solution of both systems is unique and given by the Peano-Baker series. But this series is of a little practical use. Only in the case of certain classes of LTV-system the Peano-Baker series reduces to a simpler applicable form and the solution can be found. These classes are, for example, a commutative LTV-system, LTV-system with a diagonal dynamic matrix and a periodic LTV-system, see
j=0
is in general the function of the system matrices b (t) and A (t) as well as their corresponding time derivatives. Assuming ∂Y Φ1 (yh ) = 0, (9) ∂Y2 Φ1 (yh ) = 0, LTV-transformation (8) is reduced to LTI-transformation. In other words, all the corresponding time derivatives of b (t) and A (t) disappear, and the transformation itself is
3
Proceedings of the 9th MATHMOD 34 Vienna, Austria, February 21-23, 2018
A. Galkina et al. / IFAC PapersOnLine 51-2 (2018) 31–36
with the step number k = 0, ..., N and the corresponding system dimensions nx = 4 and ny = 2. The dynamic of the resulting subsystems ˆ1 h ˆ1 hk+1 0 1 0 0 0 k ˆ1 0 1 ˆ1 0 0 1 0 h h k+1 ˆ k+2 = h ˆ 1 + 0 hk+4 ˆ1 0 0 0 1h k+2 k+3 0 0 0 0 1 ˆ1 ˆ1 h h k+3 k+4 yˆh,k+1 0 1 yˆh,k 0 . = + yˆ 0 0 1 h,k+2 yˆh,k+2 yˆh,k+1 ˆ 1 and y ˆ k . Hence, the is the simple shift in the sequences h k extended system states T ˜ k = xc,k q¯k1 x˙ c,k q¯˙k1 Fx,k , x
e.g Zadeh and Desoer (1963), Shiobara and Hori (2008). The considered LTV-subsystem of (11) with substitution of a linear function of t for yh,d is neither a commutative one, no dsolve of Maple 2017 can provide the solution. Hence the solution is not available. In the following, an alternative approximate discretization method is proposed, namely the exact time discretization of a flat space. In accordance with this method, LTV-system is transformed into a Brunovsky canonical form and then the integrator chain, that is LTI-system, is discretized exactly. Of course, the exact time discretization of the LTI-system does not correspond in general to the exact time discretization of the original LTV-system. But for the sufficiently small Ta this simplification is a reasonable one.
ˇ k = [ yh,k y˙ h,k Fy,k ] y are parametrized by the flat output and we have ˜1, ˜ k = TX (yh,d,k ) h x k ˜k, ˇ k = TY y y where 1 1 T ˜1 = h ˆ1 ˆ1 ˆ ˆ1 ˆ h h , h h h k T
4.2 System Discretization and Parameterization by Flat Output To transform the LTV-subsystem of the simplified system (11) into Brunovsky canonical form, firstly, the transformation into controllability canonical form is carried out by ¯ 1 = Tx,1 (yh ) x ¯ , which results in h ˙h ¯1 = A ¯1 + b u (y ) h ¯ Rx,1
h,d
k
R,2
T with the new state h = h h h h1,(3) , where h1 is the flat output of the LTV subsystem of the simplified system (11). Secondly, we transform the system input by the relation nx u ¯ = i=0 ai (yh,d ) h1,(i) with the coefficients of the characteristic polynomial ai (yh,d ) of ARx,1 (yh,d ) into the Brunovsky canonical form and get an integrator chain 1 h1 h˙ 0 0100 h ¨ 1 0 0 1 0 h˙ 1 0 1,(4) 1 + h . (12) ¨ h1,(3) = 0 0 0 1 h 0 1 0000 h1,(3) h1,(4) ¯1
1
˙1
¨1
k+1
k
k+2
k+3
k+4
and the input sequence u[0,N −1] = [ u0 ... uN −1 ] . Hence, the state zN depends on N entries of u[0,N −1] . Comparing with the described above parameterization, it can ˜ k and y ˜ k are be seen, that the extended system states x parametrized by only 5 and 3 entries of the flat outputs, respectively. Due to such parameterization the Jacobian matrix of the corresponding optimization problem has a sparse band structure, for which efficient linear solvers are available. The implementation of this parameterization can be found in Galkina et al. (2016) and Galkina et al. (2017). To accelerate the online calculations the transformation sequence TX (yh,d,k ) is calculated offline. T
Recall, that the LTI-subsystem of (11) is already in Brunovsky canonical form. Then we discretize both subsystems by the exact discretization with the sampling time Ta and get 4 Ta T2 T3 1 Ta 2a 6a 24 Ta3 Ta2 1 1 ¯ ¯ 0 1 Ta 2 h 6 h1 h , k+1 = k + Ta2 k+4 0 0 1 Ta 2 0 0 0 1 Ta T2 a 1 Ta ¯ k+1 = ¯ k + 2 yh,k+2 . y y 0 1 Ta Finally, the discretized subsystems are transformed into Brunovsky canonical form one more time, using the transformation into the controllability canonical form, namely ¯1, ˆ 1 = Tx,2 h h k k (13) ¯k, ˆ k = Ty,2 y y resulting in ˆ 1 = ARx,2 h ˆ 1 + bRx,2 h1 h
5. MODEL PREDICTIVE CONTROL In this section, we introduce the model predictive control (see for example, Mayne et al. (2000)) for the stabilization of the offline generated trajectory. The crucial point for the MPC implementation is its computational efficiency. It is shown in Galkina et al. (2017) that MPC in the form of a linear program (LP) can be much faster comparing to a quadratic program (QP). Moreover, the equivalent formulation of LP and QP results in the same states dynamic. In the literature, see e.g. Agrawal et al. (2017), the following definition of equivalent convex optimization problems can be found: two problems are equivalent if a solution for one can be readily converted to a solution for the other.
k+4
ˆ k+1 = ARy,2 y ˆ k + bRy,2 yh,k+2 y and the transformation of the input by the relations nx ˆ1 h1k+4 = i=0 ai (ARx,2 ) h k+i (14) ny yh,k+2 = i=0 ai (ARx,2 ) yˆh,k+i
k+1
˜ k = [ yˆh,k yˆh,k+1 yˆh,k+2 ] . y To explain the motivation behind the second transformation set (13) and (14), let us consider the general discretetime LTI-system zk+1 = Azk + buk . The solution of such system for some k = N can be written in the following form zN = AN z0 + HN u[0,N −1] with the time-varying Hankel matrix HN = AN −1 b ... Ab b
4
Proceedings of the 9th MATHMOD Vienna, Austria, February 21-23, 2018
A. Galkina et al. / IFAC PapersOnLine 51-2 (2018) 31–36
35
values corresponding to the optimal trajectory at the first step. Then the optimization vector can be initialized with the optimal solution of the previous step. MPC strategy enables the simple incorporation of the dead-time due to system discretization. To avoid the control input delay of 20 ms, the next following sample Fx (ˆ p1 , 2) and Fy (ˆ p2 , 2) of the corresponding optimal control sequences are applied.
We restrict ourselves on QP and LP. By the equivalent formulation of LP and QP we mean usage of the same constraints in an optimization program and the same structure of the cost function, but with different norms. The optimal solutions for both formulations are equal. For example, the cost function of LP in the following form JLP = δ1 + ... + δm with the positive variables δ1 , ..., δm and m ∈ N corresponds to the equivalent cost 2 function JQP = δ12 + ... + δm of QP. Different classes of convex problems fit into a spacial hierarchy, see e.g. Agrawal et al. (2017). In accordance with this hierarchy, any linear program can be considered as more general class of quadratic program. From LP to QP problem classes, both generality and difficulty increases and it is preferable to solve linear program problems using linear solvers and not quadratic ones, although in linear MPC, cost function is often a quadratic one. An optimal solution of LP is always located at ’corner points’ (intersection of constraints) of a constraint polytope. In contrast, an optimal solution of QP can be located anywhere inside the constraint polytope. Bearing the above discussion in mind, we prefer LP for the MPC design.
6. SIMULATION RESULTS To verify the implementation of the designed MPC, the numerical simulation with Matlab 2015a Simulink is presented. We consider the robustness of MPC to parameter uncertainties and external disturbances. For this aim, most of the system parameters get the deviation of ∼10%. To include the external disturbances, the slightly varying uncompensated part of the Coulomb friction force, which is independent of the velocity, is added to the optimal driving force Fx and a white noise signal is added to the carriage velocity x˙ c . The simulation results for the nonlinear approximated system (5) with the feedforward control are presented on the left side in Fig.3. It can be seen, that the payload position has the considerable error. In contrast, the system with MPC, which is realized by the open-source optimal solver lp solve for Matlab and have the prediction horizon N = 40, demonstrates almost no error. At the same time the computational time, that is, the time needed for solving the optimization tasks at each step, is lower than the sampling time Ta = 20 ms. The open-source solver lp solve implements the revised Simplex method, a scheme for ordering the computations required for the Simplex method, so that unnecessary calculations are avoided. Therefore, this method is a better choice for solving the presented problem (see, for example, Luenberger and Ye (2008)). The choice of the |∆xmax |, |∆Mbmax |, |∆Fxmax | and |∆yhmax |, |∆y˙ hmax | c depends significantly on the uncompensated friction force. The higher is the friction, the higher acting force deviation |∆Fxmax | must be set up. We also compared the performance of the open-source solver lp solve and optimal Matlab 2015a solvers, namely linprog and quadprog, for LP (15), (16) and the equivalent QP, respectively, see Table 1. The implementation of ’simplex’ or ’dual-simplex’ algorithm with linprog does not provide an appropriate solution at all.
The optimization task is separated into two independent optimization programs, for LTV- and LTI-subsystems of (11), respectively. To minimize the absolute values of the predicted errors, the new variables 1 , 2 , 3 and ε1 , ε1 are introduced. These variables are positive and restricted by the corresponding absolute values of the chosen maximal allowed errors |∆xmax |, |∆Mbmax |, |∆Fxmax | and c max max |∆yh |, |∆y˙ h |, respectively. The considered system has the finite run time horizon NT and it must be taking into account. At first, we consider the case, when NT > N with the prediction horizon N and formulate the optimization task in the form of two independent linear programs, namely minpˆ1 1 + 2 + 3 subject to p1 , k) − xc,d (k) | ≤ 1 k = 1, .., N |xc (ˆ |Mb (ˆ p1 , k) − Mb,d (k) | ≤ 2 k = 1, .., N |Fx (ˆ p1 , k) − Fx,d (k) | ≤ 3 k = 1, .., N (15) 0 ≤ 1 ≤ |∆xmax | c max 0 ≤ 2 ≤ |∆Mb | 0 ≤ 3 ≤ |∆Fxmax | and minpˆ2 ε1 + ε2 subject to p2 , k) − yh,d (k) | ≤ ε1 k = 1, .., N |yh (ˆ |y˙ h (ˆ p2 , k) − y˙ h,d (k) | ≤ ε2 k = 1, .., N (16) 0 ≤ ε1 ≤ |∆yhmax | max 0 ≤ ε2 ≤ |∆y˙ h | with the optimization vectors T ˆ1 ˆ1 . . . h pˆ1 = 1 2 3 h , 5 N −1 T pˆ2 = [ ε1 ε2 yˆh,3 . . . yˆh,N −1 ] and planned trajectories (xc,d , Mb,d , Fx,d ), (yh,d , y˙ h,d ). ˆ 1 and y ˆ 0 depend on the initial state of The states h 0 SMC and can be calculated offline. When NT = N the prediction horizon N must be replaced by the decreasing prediction horizon N − n with the new actual step n = 0...N − 1. In general the initial value for the optimization vectors can be set zero. To achieve higher calculation efficiency, the optimization vectors are initialized by the
Optimization Tools quadprog linprog lp solve
Algorithm active-set active-set revised simplex
Calc. Time (s) ∼0.12 ∼0.06 <0.02
Table 1: Comparison of calculation time of optimal solvers implemented in Matlab 2015a CONCLUSION In this contribution, the tracking problem for the flexible SMC was solved by means of MPC. It was shown, that the simplified linear model of SMC can be used, to design an efficient model predictive control in the form of a linear program applying flatness-based parameterization. The parameterization of the system by the special flat output 5
Proceedings of the 9th MATHMOD 36 Vienna, Austria, February 21-23, 2018
XY - payload position
0.2 0
0.2
0.4 0.6 X (m) Mb - bending moment
Mb (Nm)
1
Trajectory Feedforward
0
0.5
1
1.5
2
20 0 -20
0.2
0.4
0.6 X (m) Mb - bending moment
0.8
1
Trajectory MPC
0
0.5
60 Trajectory Feedforward
40
0
0
-1
2.5
t (s) Fx - driving force
60
0.2
1
0
-1
Trajectory MPC
0.4
0
0.8
Mb (Nm)
0
Y (m)
Trajectory Feedforward
0.4
XY - payload position
0.6
Fx (N)
Y (m)
0.6
Fx (N)
A. Galkina et al. / IFAC PapersOnLine 51-2 (2018) 31–36
1
1.5 t (s) Fx - driving force
2
2.5
Trajectory MPC
40 20 0
0
0.5
1
1.5
2
2.5
t (s)
-20
0
0.5
1
1.5
2
2.5
t (s)
Fig. 3. Simulation results: left - approximate nonlinear model with feedforward control; right - approximate nonlinear model with MPC facilitates the formulation of an optimization task and provides the band structure of the corresponding Jacobian matrix. Such optimization task can be simply solved by means of any open-source linear optimal solver. The simulation of the designed MPC together with the nonlinear approximated model of SMC subjected to model uncertainties and external disturbances, such as uncompensated friction force, demonstrated a robust performance. The prediction horizon of MPC can be chosen high enough to avoid stability problems. For the future work, the implementation of MPC with the laboratory model of SMC and is planned. Also, the design of the optimal trajectory using the presented flatness-based parameterization will be investigated in the future.
Galkina, A., Gafur, I., and Schlacher, K. (2017). Model predictive control with linear programming for the strip infeed in hot rolling mills. In Proceedings of the 20th World Congress of the International Federation of Automatic Control, 1–6. Hajdu, S. and Gaspar, P. (2014). From modeling to robust control design of single-mast stacker cranes. Acta Polytechnica Hungarica, 11(10), 135–149. Luenberger, D. and Ye, Y. (2008). Linear and Nonlinear Programming. Springer, 233 Spring Street, New York, NY 10013, USA, 4 edition. Mayne, D., Rawlings, J., Rao, C., and Scokaert, P. (2000). Constrained model predictive control: Stability and optimality. Automatica, 36, 789–814. Rams, H., Sch¨oberl, M., and Schlacher, K. (2017). Optimal control of a single mast stacker crane. IEEE Transactions on Control System Technology. Shiobara, H. and Hori, N. (2008). Numerical exact discrete-time-model of linear time-varying systems. In Proceedings of International Conference on Control, Automation and Systems, 2314–2318. Staudecker, M., Schlacher, K., and Hansl, R. (2008). Passivity based control and time optimal trajectory planning of a single mast stacker crane. In Proceedings of the 17th IFAC World Congress. Zadeh, L. and Desoer, C. (1963). Linear System Theory: The State Space Approach. McGraw-Hill, New York.
REFERENCES Agrawal, A., Verschueren, R., Diamond, S., and Boyd, S. (2017). A rewriting system for convex optimization problems. https://arxiv.org/abs/1709.04494. Aschemann, H. and Schindele, D. (2010). Model predictive trajectory control for high-speed rack feeders. In In Model Predictive Control, 183–198. InTech. Bachmayer, M., Ulbrich, H., and Rudolph, J. (2011). Flatness-based control of a horizontally moving erected beam with a point mass. Mathematical and Computer Modelling of Dynamical Systems, 17(1), 49–69. Galkina, A., Gafur, I., and Schlacher, K. (2016). Modellpr¨ adiktive Regelung f¨ ur das Bandeinf¨ adeln in der Warmwalzstrasse. BHM Berg- und Huettenmaennische Monatshefte, 161(11), 520—-525. doi:10.1007/s00501016-0539-6. 6
ScienceDirect
IFAC PapersOnLine 51-2 (2018) 31–36
Flatness-Based Model Predictive Control Flatness-Based Model Predictive Control Flatness-Based Model Predictive Control Flatness-Based Model Predictive Control with Linear Programming for a Single Mast with Linear Programming for a Single Mast with Linear Programming for a Single Mast with Linear Programming for a Single Mast Stacker Crane Stacker Crane Stacker Crane Stacker Crane ∗ ∗
Galkina, A. ∗ Schlacher, K. ∗ Galkina, A. ∗∗ Schlacher, K. ∗∗ Galkina, Galkina, A. A. Schlacher, Schlacher, K. K. ∗ ∗ Institute of Automatic Control and Control Systems Technology, ∗ Institute of Automatic Control and Control Systems Technology, of Control Control Systems Technology, Johannes University Linz, and Altenberger 4040 Linz, ∗ Institute InstituteKepler of Automatic Automatic Control and Control Strasse Systems69, Technology, Johannes Kepler University Linz, Altenberger Strasse 69, 4040 Linz, Johannes Kepler University Linz, Altenberger Strasse 69, Austria,(e-mail: [email protected]). Johannes Kepler University Linz, Altenberger Strasse 69, 4040 4040 Linz, Linz, Austria,(e-mail: [email protected]). Austria,(e-mail: [email protected]). Austria,(e-mail: [email protected]). Abstract: This paper deals with a trajectory tracking model predictive control for a single Abstract: paper deals with a trajectory model predictive control for aa single Abstract: This paper deals with a tracking model predictive control for mast stackerThis crane, which is used automatic tracking storage or retrieval of payloads automated Abstract: This paper deals with for a trajectory trajectory tracking model predictive controlin for a single single mast stacker crane, which is used for automatic storage or retrieval of payloads in automated mast stacker crane, which is used for automatic storage or retrieval of payloads in automated warehouses. The mathematical model of the plant is a distributed parameter one, but it admits mast stackerThe crane, which is used for of automatic storage or retrieval of payloads in automated warehouses. mathematical model the plant is a distributed parameter one, but it admits warehouses. mathematical the is one, it an excellent The approximation by model a flat of lumped parameter system. parameter The approximate nonlinear warehouses. The mathematical model of the plant plant is aa distributed distributed parameter one, but but nonlinear it admits admits an excellent approximation by a flat lumped parameter system. The approximate an excellent approximation by a flat lumped parameter system. The approximate nonlinear system is simplified further and the resulting linear time-varying system is implemented to an excellent approximation by a flat lumped parameter system. The approximate nonlinear system is simplified further and the resulting linear time-varying system is implemented to system is simplified further and the resulting linear time-varying system is implemented to design a model predictive controller. To provide an efficient formulation of an optimization task system ismodel simplified further and theToresulting linear time-varying system implementedtask to design predictive controller. provide an efficient formulation of anisoptimization design a model controller. To provide formulation of task for thea the LTV-model parametrized by a special flat output. Then, the model design acontroller, model predictive predictive controller.is To provide an an efficient efficient formulation of an an optimization optimization task for the controller, the LTV-model is parametrized by a special flat output. Then, the model for the LTV-model is by special output. Then, the predictive control isthe formulated in the form of a linear For solving derived linear for the controller, controller, the LTV-model is parametrized parametrized by aaprogram. special flat flat output. the Then, the model model predictive control is formulated in the form of aa linear program. For solving the derived linear predictive control is formulated in the form of linear program. For solving the derived linear program, an open-source optimal solver lp solve is chosen. To reduce a mechanical stress at predictive control is formulated in solver the form of a linear program. For solving the derived linear program, an open-source optimal lp solve is chosen. To reduce a mechanical stress at program, an open-source optimal solver lp solve is chosen. To reduce a mechanical stress clamping of a mast and an overall accessibility time, a special optimal trajectory calculated program, an open-source optimal solver lp solvetime, is chosen. To reduce a trajectory mechanicalcalculated stress at at clamping of a mast and an overall accessibility aa special optimal clamping of a mast and an overall accessibility time, special optimal trajectory calculated offline is chosen for the tracking by the designed controller. Finally, the simulation results of clamping of a mast and tracking an overall accessibility time, a special optimal trajectory calculated offline is chosen for the by the designed controller. Finally, the simulation results of offline is chosen for the tracking by the designed controller. Finally, the simulation results the nonlinear approximate model, subjected to model uncertainties and external disturbances, offline is chosen for the tracking by the designed controller. Finally,and theexternal simulation results of of the nonlinear approximate model, subjected to model uncertainties disturbances, the nonlinear approximate model, subjected model uncertainties and external disturbances, with the designed model predictive controllerto are presented. the nonlinear approximate model, subjected to model uncertainties and external disturbances, with the designed model predictive controller presented. with the designed model controller are are presented. with theIFAC designed model predictive predictive are presented. © 2018, (International Federation controller of Automatic Control) Hosting by Elsevier Ltd. All rights reserved. Keywords: Model predictive control, Flatness, Linear programming, Flexible construction Keywords: Model predictive control, Flatness, Linear programming, Flexible construction Keywords: Keywords: Model Model predictive predictive control, control, Flatness, Flatness, Linear Linear programming, programming, Flexible Flexible construction construction 1. INTRODUCTION of modeling, see e.g. Bachmayer et al. (2011), considers the 1. INTRODUCTION INTRODUCTION of modeling, modeling, see e.g.during Bachmayer et al. al. (2011), considersThe the 1. of e.g. Bachmayer et considers the payload beingsee fixed positioning of the carriage. 1. INTRODUCTION of modeling, see e.g.during Bachmayer et al. (2011), (2011), considersThe the payload being fixed positioning of the carriage. being during positioning of the carriage. The second type of fixed the modeling approach assumes the payload A flexible mechanical system, depicted in Fig.1, represents payload payload being fixed during positioning of the carriage. The secondmoved type of of the the its modeling approach assumes the payload payload flexiblemast mechanical system, depicted in type modeling assumes the being positionapproach is treated as an uncertainty. A flexible mechanical depictedalso in Fig.1, Fig.1, represents aA single stacker system, crane (SMC), called represents a flexible second second type ofand the its modeling approach assumes the payload being moved and position is treated as an uncertainty. A flexiblemast mechanical system, depicted in Fig.1, represents a single stacker crane (SMC), also called a flexible being moved and its position is treated as an uncertainty. In Hajdu andand Gaspar (2014) this approach with arack single mast stacker crane also called flexible feeder. is usually used for automatic being moved its position is treated as antogether uncertainty. Hajdu and modeling Gaspar (2014) this approach together with a single mastSMC stacker crane (SMC), (SMC), also called aa storage flexible In rack feeder. SMC is usually usually used for for automatic storage Hajdu and Gaspar (2014) this approach together with aIn multi-body technique is presented. The third rack feeder. SMC is used automatic storage or retrieval of payloads in automated warehouses. The In Hajdu and Gaspar (2014) this approach together with a multi-body modeling technique is presented. The third rack feeder. SMC is usually used for automatic storage or retrieval retrieval of payloads payloads in initial automated warehouses. The atype multi-body modeling technique is presented. The third third of the modeling approach exploiting a fast underlying or of in automated warehouses. The payload is moved from the position to the desired a multi-body modeling technique is presented. The type of of the the modeling approach exploiting fast underlying or retrieval of payloads in initial automated warehouses. The type payload is moved moved from the the position toreduced the desired desired modeling exploiting aaa fast underlying control on approach the current converter both for the payload is from position to the end position. The weight of initial the mast is often due velocity type of the modeling approach exploiting fast underlying velocity control on the current converter both for and the payload is moved from the initial position toreduced the desired end position. The weight of the mast is often due velocity control on the current converter both for the payload and carriage cancurrent be found in Aschemann end position. The weight of the mast is often reduced due to economic reasons. Therefore, vibrations are generated velocity control on the converter both for and the payload and carriage can be found in Aschemann end position. The weight of the mast is often reduced due to economic economic reasons. Therefore, Therefore, vibrations areHence, generated payload and carriage can be found in Aschemann and Schindele (2010). Another modeling approach is presented to reasons. vibrations are generated during transportation and positioning process. the payload and carriage can be found in Aschemann and Schindele (2010). Another modeling approachassumes is presented presented to economic reasons. Therefore, vibrations areHence, generated during transportation and positioning positioning process. (2010). modeling approach is in Staudecker et Another al. (2008). This approach the during transportation and process. Hence, the Schindele structural vibrations must be taking into account in the Schindele (2010). Another modeling approachassumes is presented in Staudecker et al. (2008). This approach the during transportation and positioning process. Hence, the structural vibrations must be taking into account in the in Staudecker et al. (2008). This approach assumes the controlled motion of the payload during the carriage transstructural vibrations must be taking into account in the model being derived for control of the positioning process. in Staudecker et al. (2008). This approach assumes the controlled motion of the the payload during the carriage carriage transstructural vibrations must be taking into accountprocess. in the controlled model being derived for control of the positioning motion of payload during the transportation. The first and second approaches do not allow model control of the positioning process. SMC isbeing an derived essentialfor element of automated warehouses controlled motion of the payload during the carriage transportation. The first and second approaches do not allow model being derived for control of the positioning process. SMCcontributes is an an essential element of automated automated warehouses The first second not minimizing the overall accessibility time of thedo storage and SMC is of warehouses and a lotelement to its overall accessibility time. portation. portation. The first and and second approaches approaches not allow allow minimizing the overall overall accessibility time of of the thedostorage storage and SMC is an essential essential element of automated warehouses and contributes a lot to its overall accessibility time. minimizing the accessibility time and retrieval processes. and contributes aa lot to overall accessibility time. Therefore, the control to position the system as minimizing the overall accessibility time of the storage and retrieval processes. and contributes lotpurpose to its its is overall accessibility time. Therefore, the control purpose is to position the system as processes. Therefore, the purpose is position system as retrieval fast as possible. To achieve this purpose thethe feedforward retrieval processes. This paper deals with MPC for stabilization of the offline Therefore, the control control purpose is to to position system as fast as as possible. possible. To achieve achieve this purpose thethe feedforward This paper deals with with MPC for for stabilization ofthe themodel, offline fast To this purpose the feedforward control design based on modern control concepts such as This paper deals MPC stabilization the offline generated optimal trajectory and is based onof fast as possible. To achieve this purpose the feedforward control design design optimal based on ontrajectory modern control control concepts such as as This paperoptimal deals with MPC for stabilization ofthe themodel, offline generated trajectory and is based on control based modern concepts such flatness-based planning can be applied. optimal trajectory and on the model, which is presented in Staudecker etis al.based (2008). The mathecontrol design optimal based ontrajectory modern control concepts such as generated flatness-based planning can be applied. generated optimal trajectory and is based on the model, which is presented in Staudecker et al. (2008). The matheflatness-based optimal trajectory be applied. applied. But the considered system iscan subjected to maticalis model presented in Staudecker et The of the is a distributed parameter one, flatness-based optimalmechanical trajectory planning planning can be But the considered considered mechanical system is to which which is model presented in plant Staudecker et al. al. (2008). (2008). The mathemathematical of excellent the plant is aa distributed distributed parameter one, But the mechanical system is subjected subjected to parameter uncertainties and external disturbances such as matical model of the plant is parameter one, but it admits an approximation by a flat lumped But the considered mechanical system is subjected to parameter uncertainties and external disturbances such as as matical modelan of excellent the plantapproximation is a distributedbyparameter one, but it admits a flat lumped parameter uncertainties and external disturbances such friction force. Therefore, feedback extension is desired, to but it admits an excellent approximation by aa flat lumped parameter system. To reduce the mechanical stress due parameter uncertainties and external disturbances such as friction force. Therefore, feedback extension is desired, to but it admits an excellent approximation by flat lumped parameter system. To reduce the mechanical stress due friction force. Therefore, feedback extension is to stabilize the planned offline trajectory. Stabilization of the system. To reduce the mechanical stress due to a bending moment acting on the mast clamping and friction force. Therefore, feedback extension is desired, desired, to parameter stabilize the planned offline trajectory. Stabilization of the parameter system. To reduce the mechanical stress due to a bending moment acting on the mast clamping and stabilize the offline Stabilization of desired trajectory is the maintrajectory. objective of the contribution. aareduce bending acting on clamping and the moment transporting the mast optimal trajectory stabilize the planned planned offline Stabilization of the the to desired trajectory isthe the maintrajectory. objective of the contribution. contribution. to bending moment actingtime, on the the mast clamping and to reduce the transporting time, the optimal trajectory desired trajectory is the main objective of the For this purpose, design of a flatness-based model to reduce the transporting time, optimal described et al. (2017) is chosen to trajectory be stabidesired isthe the design main objective of the contribution. For this thistrajectory purpose, of linear flatness-based model to reduce by the Rams transporting time, the the optimal trajectory described by Rams etminimizes al. (2017) (2017) is bending chosen tomoment be stabistabiFor purpose, design of aa flatness-based model predictive control the (MPC) with programming is described by Rams et al. is chosen to be lized. This trajectory the at For this purpose, the design of a flatness-based model predictive control control (MPC) (MPC) with with linear linear programming programming is is described by Rams etminimizes al. (2017)the is bending chosen tomoment be stabilized. This trajectory at predictive proposed. lized. This trajectory minimizes the bending moment at clamping of the mast and, at the same time, moves predictive control (MPC) with linear programming is the proposed. lized. This trajectory minimizes the bending moment at the clamping of the mast and, at the same time, moves proposed. clamping of mast and, the time, moves payload from starting position to the poproposed. There may be found several approaches both for mathe- the the clamping of the thethe mast and, at at the same same time,end moves the payload from the starting position to the end poThere may may be found found several approaches both both for mathepayload the to position as fastfrom as possible. It is position worth mentioning, There be approaches mathematical modeling and several for feedforward/feedback control of the the payload the starting starting to the the end end that position as fast fastfrom as this possible. It is position worth mentioning, that There may be found several approaches both for for mathematical modeling and for feedforward/feedback control of sition as as possible. It is worth mentioning, that the presented in paper MPC concept can be applied matical modeling and for feedforward/feedback control of SMC in the literature. The modeling approaches differ in sition as fast as possible. It is worth mentioning, that the presented in this paper MPC concept can be applied matical and for feedforward/feedback control of the SMC in modeling the payload literature. The modeling approaches differ in presented in paper MPC concept can be applied with any trajectory could be suitable SMC in the literature. approaches differ the way andThe the modeling mast are treated, as well as in the in this thisand paper MPC concept for canthe be controlapplied withpresented any trajectory trajectory and could be suitable suitable for the controlSMC in the the payload literature. The modeling approaches differ in with the way the and the mast are treated, as well as in any and could be for the controlthe way the payload and the mast are treated, as well as in with any trajectory and could be suitable for the controlphysical variables being used as control inputs. One type the way the payload and the mast are treated, as well as in physical variables variables being used used as control control inputs. One One type physical physical variables being being used as as control inputs. inputs. One type type
Copyright © 2018, 2018 IFAC 1 Hosting by Elsevier Ltd. All rights reserved. 2405-8963 © IFAC (International Federation of Automatic Control) Copyright © 2018 IFAC 1 Copyright ©under 2018 responsibility IFAC 1 Control. Peer review of International Federation of Automatic Copyright © 2018 IFAC 1 10.1016/j.ifacol.2018.03.006
Proceedings of the 9th MATHMOD 32 Vienna, Austria, February 21-23, 2018
xk
A. Galkina et al. / IFAC PapersOnLine 51-2 (2018) 31–36
point mass mh moves along the beam centerline. The forces Fx and Fy acting on the carriage and the lifting unit, respectively, serve as control inputs.
mk
w(Y, t) xh
yh
Concerning a mathematical notation, the following shortcuts are used, namely ∂t = ∂/∂t and ∂Y = ∂/∂Y as well as ∂tk = ∂ k /∂tk and ∂Yk = ∂ k /∂Y k , k = 2, ..., n. When dealing with functions u = u(t), it is written u˙ = ∂t u, u ¨ = ∂t2 u and u(i) for the time derivatives of the higher order.
mh
Fy
Y X
Fx
2.1 System Modeling
mc
To shorten some expressions, we interchange sometimes the function, e.g. w (Y, t), and the coordinate, e.g. w, notation. Taking the mentioned above assumptions into account, the kinetic energy of the system is ˆL mk 2 ρA mc 2 mh 2 2 2 Ek = x˙ + x˙ h + y˙ h + x˙ + (∂t w) dY 2 c 2 2 k 2
xc Fig. 1. Single mast stacker crane oriented model presented, for example, in Aschemann and Schindele (2010) as well. The paper consists of six main parts. At first, to keep a consistency of the proposed approach, modeling of the infinite-dimensional system and the system approximation by a finite-dimensional model is recapitulated in Section 2. In general, the nonlinear approximate model can be used together with a nonlinear MPC for stabilization of a desired trajectory. But the optimization problem will be a nonlinear one, that increases the calculation time considerably and it is hardly applicable to a real-time system. For this reason, in Section 3, the flatness analysis, which is proposed in Staudecker et al. (2008), is implemented to simplify the approximate model. As a result, a linear timevarying (LTV) model is derived. We consider the problem of stabilization of the optimal trajectory being generated offline and solve it by means of model predictive control. For MPC design the simplified LTV-model must be discretized with respect to time. For this aim, discretization of a flat space is proposed in Section 4. Also, to improve an efficiency of MPC realized by the free open-source optimal solver lp solve, in Section 4 the discretized flat system states are parameterized by a flat output again. Afterward, in Section 5, the stabilizing MPC in the form of a linear program is presented. The flatness-based MPC design with linear programming exploiting the simplified system model of SMC represents the main contribution of this paper. Finally, the simulation results with the nonlinear finite-dimensional system, which is subjected to parameter uncertainties and external disturbances, are presented in Section 6.
0
and the potential energy is
1 Ep = mh gyh + EI 2
ˆL 0
∂Y2 w
2
dY.
The contribution of the external forces can be expressed as Eext = xc Fx + yh Fy . Hence, the Lagrangian of the system can be written as L = Ekin − Epot + Eext . To derive the equations of motion for the SMC the Hamilton’s principle is applied and the following equations of motion are derived in the form of one partial differential equation ρA∂t2 w = −EI∂Y4 w
and the set of ordinary differential equations mc x ¨c = Fx − EI∂Y3 w (0, t) mh y¨h = Fy − ¨h ∂Y w (y h , t) mh g − mh x mh x ¨h = EI ∂Y3 w yh− , t − ∂Y3 w yh+ , t ¨k = EI∂Y3 w (L, t) . mk x Also, the boundary conditions must be satisfied w (0, t) = xc , w (L, t) = xk , ∂Y w (0, t) = 0 , ∂Y2 w (L, t) = 0 , w (yh , t) = xh , ∂Y2 w yh− , t = ∂Y2 w yh+ , t .
(1)
(2)
(3)
For modeling details see Staudecker et al. (2008). The values yh+ and yh− denote the most upper and the most lower position yh of the payload. Next, the mixed-dimensional system equations (1), (2) and (3) are discretized in accordance with the Rayleigh-Ritz method and the finitedimensional approximation is derived.
2. MODELING In this section, we recapitulate briefly the derivation of the mathematical model, which is proposed in Staudecker et al. (2008). The carriage, the lifting unit and the tip with the corresponding masses mc , mh and mk are modeled as rigid bodies with the mechanical degrees of freedom xc , xh , yh and xk in the fixed coordinate frame XY (see Fig. 1). The displacement of the mast centerline is denoted by w (Y, t) and it is assumed that the mast can be considered as a beam, which satisfies the Euler-Bernoulli hypothesis, having length L, uniform mass density ρ, uniform cross section surface A and uniform flexural rigidity EI. The
2.2 System Approximation Applying the Rayleigh-Ritz method, a finite-dimensional approximation of the system model, represented by (1), (2) and (3), can be derived. According to this method the mast deflection w (Y, t) is approximated by the first-order Ritz ansatz function w∗ (Y, t) = xc (t) + Φ1 (Y ) q¯1 (t) (4) 2
Proceedings of the 9th MATHMOD Vienna, Austria, February 21-23, 2018
A. Galkina et al. / IFAC PapersOnLine 51-2 (2018) 31–36
only the function of yh,d . Substituting the simplification (9) into (7) yields the LTI-system 0 1 0 ˙y ¯= ¯+ y v¯ (10) 0 0 1
with the new generalized coordinate q¯1 , and the appropriately chosen spatial basis function 4
with the input v¯ = (Fy − gmh )/mh . Hence, the resulting simplified system can be written in the following form ¯˙ = Ax (yh,d ) x ¯ + bx (yh,d ) u x ¯ (11) ¯ + by v¯ ¯˙ = Ay y y with the LTV- and LTI-subsystems. The numerical verification (see Fig.(2)) shows, that applying the same input to the approximate (5) and simplified approximate model (11), the only small difference in the state dynamics occurs. Hence, the considered simplified model (11) can be used for the design of MPC in the form of a linear program. XY - payload position
0.6
Y (m)
Then, based on the flatness analysis, which is also described in Staudecker et al. (2008), the simplification of the approximate system (5) is carried out.
0.4 0.2 0
3. SYSTEM SIMPLIFICATION
0
0.2
0.4
0.6
0.8
Mb - bending moment
1
Exact AM Simplified AM
Mb (Nm)
3
Φ1 (Y ) = 6 (Y /L) − 4 (Y /L) + (Y /L) . Then, using the Ritz ansatz function (4) and applying the Hamilton’s principle, the approximate equations of motion can be written in the form of the explicit nonlinear ordinary differential equation of the second order ¨ + C (q, q) ˙ = Gu. M (q) q (5) The derived approximate system has the generalized co T ordinate q = xc q¯1 yh and the control input u = T [ Fx Fy ] . The bending moment acting at the clamping of the mast is proportional to q¯1 and can be expressed as Mb = 12EI q¯1 /L2 . The details and verification of the considered approximation can be found in Staudecker et al. (2008).
0.5 0 -0.5
1
-1
Exact AM Simplified AM
0
X (m) y˙ h - payload velocity
0.4 y˙ h (m/s)
In spite of the fact, that the state space representation of (5) is not state feedback equivalent to a linear timeinvariant (LTI) system, nonetheless it is flat. In Staudecker et al. (2008) it is shown, how the system decomposition can be used for the flatness analysis of the approximate system (5). It is also shown that yh is a flat output corresponding to the vertical dynamics and the desired trajectory yh,d (t) with the corresponding time derivatives are known (are computed offline). The system (5) is decomposed into a linear time-varying and nonlinear subsystems, namely ¯˙ = A (yh,d , y˙ h,d , y¨h,d ) x ¯ + b (yh,d ) u x ¯ (6) 1 1 T ˙ ¯ = xc q¯ x˙ c q¯ with the state vector x , the input u ¯ = Fx and ¯˙ = f (¯ y y, u ˆ) (7)
2
0.2 0 -0.2 -0.4
Exact AM Simplified AM
0
0.5
1
1.5 t (s)
2
˙ b (Nm/s) M
2
33
2.5
0.5 1 1.5 2 t (s) ˙ b - bending moment rate of change M
2.5
0 -2 -4
Exact AM Simplified AM
0
0.5
1
1.5
2
2.5
t (s)
Fig. 2. Approximate model (AM) (5) and simplified approximate model (11) 4. FLATNESS-BASED PARAMETRIZATION For the formulation of an optimization problem a discretetime system representation of (11) must be found. A dynamical system can be discretized by means of an exact or non-exact discretization method. To discretize the LTVsystem and to derive the efficient formulation of an optimization problem in the form of a linear program, the approach of a flatness-based parameterization of the system states, presented in Galkina et al. (2017) is used. We transform the simplified system into Brunovsky canonical form. In other words, the transformation of LTV-system into LTI-system is used. Then the integrator chain is discretized with respect to time by means of the exact time discretization and the discretized system is transformed into Brunovsky canonical form one more time.
¯ = [ yh y˙ h ] and the input u with the state y ˆ = T 1 1 xc q¯ x˙ c q¯˙ Fx Fy . The system (6) can be treated as a linear time-varying one and the corresponding methods and approaches for LTV-system can be implemented. Since for a linear system the properties of flatness and controllability are equivalent, we may obtain the second flat output corresponding to the dynamics in X-direction by transforming (6) into a time-varying controllability canonical form. Now, for the considered approximate model of SMC further simplifications can be made, to simplify the formulation of the linear optimization program for the stabilizing MPC. It is worth mentioning that the LTVtransformation 4 Tj1 (t) x ¯j (8) h1 = T
4.1 Motivation For the exact time discretization, a solution of a dynamic system is needed. In the case of LTI-system, a general solution is always available. Conversely, LTV- and nonlinear systems do not have a general solution. Both LTI- and LTV-systems need for a transition matrix to get a solution. In general, the solution of both systems is unique and given by the Peano-Baker series. But this series is of a little practical use. Only in the case of certain classes of LTV-system the Peano-Baker series reduces to a simpler applicable form and the solution can be found. These classes are, for example, a commutative LTV-system, LTV-system with a diagonal dynamic matrix and a periodic LTV-system, see
j=0
is in general the function of the system matrices b (t) and A (t) as well as their corresponding time derivatives. Assuming ∂Y Φ1 (yh ) = 0, (9) ∂Y2 Φ1 (yh ) = 0, LTV-transformation (8) is reduced to LTI-transformation. In other words, all the corresponding time derivatives of b (t) and A (t) disappear, and the transformation itself is
3
Proceedings of the 9th MATHMOD 34 Vienna, Austria, February 21-23, 2018
A. Galkina et al. / IFAC PapersOnLine 51-2 (2018) 31–36
with the step number k = 0, ..., N and the corresponding system dimensions nx = 4 and ny = 2. The dynamic of the resulting subsystems ˆ1 h ˆ1 hk+1 0 1 0 0 0 k ˆ1 0 1 ˆ1 0 0 1 0 h h k+1 ˆ k+2 = h ˆ 1 + 0 hk+4 ˆ1 0 0 0 1h k+2 k+3 0 0 0 0 1 ˆ1 ˆ1 h h k+3 k+4 yˆh,k+1 0 1 yˆh,k 0 . = + yˆ 0 0 1 h,k+2 yˆh,k+2 yˆh,k+1 ˆ 1 and y ˆ k . Hence, the is the simple shift in the sequences h k extended system states T ˜ k = xc,k q¯k1 x˙ c,k q¯˙k1 Fx,k , x
e.g Zadeh and Desoer (1963), Shiobara and Hori (2008). The considered LTV-subsystem of (11) with substitution of a linear function of t for yh,d is neither a commutative one, no dsolve of Maple 2017 can provide the solution. Hence the solution is not available. In the following, an alternative approximate discretization method is proposed, namely the exact time discretization of a flat space. In accordance with this method, LTV-system is transformed into a Brunovsky canonical form and then the integrator chain, that is LTI-system, is discretized exactly. Of course, the exact time discretization of the LTI-system does not correspond in general to the exact time discretization of the original LTV-system. But for the sufficiently small Ta this simplification is a reasonable one.
ˇ k = [ yh,k y˙ h,k Fy,k ] y are parametrized by the flat output and we have ˜1, ˜ k = TX (yh,d,k ) h x k ˜k, ˇ k = TY y y where 1 1 T ˜1 = h ˆ1 ˆ1 ˆ ˆ1 ˆ h h , h h h k T
4.2 System Discretization and Parameterization by Flat Output To transform the LTV-subsystem of the simplified system (11) into Brunovsky canonical form, firstly, the transformation into controllability canonical form is carried out by ¯ 1 = Tx,1 (yh ) x ¯ , which results in h ˙h ¯1 = A ¯1 + b u (y ) h ¯ Rx,1
h,d
k
R,2
T with the new state h = h h h h1,(3) , where h1 is the flat output of the LTV subsystem of the simplified system (11). Secondly, we transform the system input by the relation nx u ¯ = i=0 ai (yh,d ) h1,(i) with the coefficients of the characteristic polynomial ai (yh,d ) of ARx,1 (yh,d ) into the Brunovsky canonical form and get an integrator chain 1 h1 h˙ 0 0100 h ¨ 1 0 0 1 0 h˙ 1 0 1,(4) 1 + h . (12) ¨ h1,(3) = 0 0 0 1 h 0 1 0000 h1,(3) h1,(4) ¯1
1
˙1
¨1
k+1
k
k+2
k+3
k+4
and the input sequence u[0,N −1] = [ u0 ... uN −1 ] . Hence, the state zN depends on N entries of u[0,N −1] . Comparing with the described above parameterization, it can ˜ k and y ˜ k are be seen, that the extended system states x parametrized by only 5 and 3 entries of the flat outputs, respectively. Due to such parameterization the Jacobian matrix of the corresponding optimization problem has a sparse band structure, for which efficient linear solvers are available. The implementation of this parameterization can be found in Galkina et al. (2016) and Galkina et al. (2017). To accelerate the online calculations the transformation sequence TX (yh,d,k ) is calculated offline. T
Recall, that the LTI-subsystem of (11) is already in Brunovsky canonical form. Then we discretize both subsystems by the exact discretization with the sampling time Ta and get 4 Ta T2 T3 1 Ta 2a 6a 24 Ta3 Ta2 1 1 ¯ ¯ 0 1 Ta 2 h 6 h1 h , k+1 = k + Ta2 k+4 0 0 1 Ta 2 0 0 0 1 Ta T2 a 1 Ta ¯ k+1 = ¯ k + 2 yh,k+2 . y y 0 1 Ta Finally, the discretized subsystems are transformed into Brunovsky canonical form one more time, using the transformation into the controllability canonical form, namely ¯1, ˆ 1 = Tx,2 h h k k (13) ¯k, ˆ k = Ty,2 y y resulting in ˆ 1 = ARx,2 h ˆ 1 + bRx,2 h1 h
5. MODEL PREDICTIVE CONTROL In this section, we introduce the model predictive control (see for example, Mayne et al. (2000)) for the stabilization of the offline generated trajectory. The crucial point for the MPC implementation is its computational efficiency. It is shown in Galkina et al. (2017) that MPC in the form of a linear program (LP) can be much faster comparing to a quadratic program (QP). Moreover, the equivalent formulation of LP and QP results in the same states dynamic. In the literature, see e.g. Agrawal et al. (2017), the following definition of equivalent convex optimization problems can be found: two problems are equivalent if a solution for one can be readily converted to a solution for the other.
k+4
ˆ k+1 = ARy,2 y ˆ k + bRy,2 yh,k+2 y and the transformation of the input by the relations nx ˆ1 h1k+4 = i=0 ai (ARx,2 ) h k+i (14) ny yh,k+2 = i=0 ai (ARx,2 ) yˆh,k+i
k+1
˜ k = [ yˆh,k yˆh,k+1 yˆh,k+2 ] . y To explain the motivation behind the second transformation set (13) and (14), let us consider the general discretetime LTI-system zk+1 = Azk + buk . The solution of such system for some k = N can be written in the following form zN = AN z0 + HN u[0,N −1] with the time-varying Hankel matrix HN = AN −1 b ... Ab b
4
Proceedings of the 9th MATHMOD Vienna, Austria, February 21-23, 2018
A. Galkina et al. / IFAC PapersOnLine 51-2 (2018) 31–36
35
values corresponding to the optimal trajectory at the first step. Then the optimization vector can be initialized with the optimal solution of the previous step. MPC strategy enables the simple incorporation of the dead-time due to system discretization. To avoid the control input delay of 20 ms, the next following sample Fx (ˆ p1 , 2) and Fy (ˆ p2 , 2) of the corresponding optimal control sequences are applied.
We restrict ourselves on QP and LP. By the equivalent formulation of LP and QP we mean usage of the same constraints in an optimization program and the same structure of the cost function, but with different norms. The optimal solutions for both formulations are equal. For example, the cost function of LP in the following form JLP = δ1 + ... + δm with the positive variables δ1 , ..., δm and m ∈ N corresponds to the equivalent cost 2 function JQP = δ12 + ... + δm of QP. Different classes of convex problems fit into a spacial hierarchy, see e.g. Agrawal et al. (2017). In accordance with this hierarchy, any linear program can be considered as more general class of quadratic program. From LP to QP problem classes, both generality and difficulty increases and it is preferable to solve linear program problems using linear solvers and not quadratic ones, although in linear MPC, cost function is often a quadratic one. An optimal solution of LP is always located at ’corner points’ (intersection of constraints) of a constraint polytope. In contrast, an optimal solution of QP can be located anywhere inside the constraint polytope. Bearing the above discussion in mind, we prefer LP for the MPC design.
6. SIMULATION RESULTS To verify the implementation of the designed MPC, the numerical simulation with Matlab 2015a Simulink is presented. We consider the robustness of MPC to parameter uncertainties and external disturbances. For this aim, most of the system parameters get the deviation of ∼10%. To include the external disturbances, the slightly varying uncompensated part of the Coulomb friction force, which is independent of the velocity, is added to the optimal driving force Fx and a white noise signal is added to the carriage velocity x˙ c . The simulation results for the nonlinear approximated system (5) with the feedforward control are presented on the left side in Fig.3. It can be seen, that the payload position has the considerable error. In contrast, the system with MPC, which is realized by the open-source optimal solver lp solve for Matlab and have the prediction horizon N = 40, demonstrates almost no error. At the same time the computational time, that is, the time needed for solving the optimization tasks at each step, is lower than the sampling time Ta = 20 ms. The open-source solver lp solve implements the revised Simplex method, a scheme for ordering the computations required for the Simplex method, so that unnecessary calculations are avoided. Therefore, this method is a better choice for solving the presented problem (see, for example, Luenberger and Ye (2008)). The choice of the |∆xmax |, |∆Mbmax |, |∆Fxmax | and |∆yhmax |, |∆y˙ hmax | c depends significantly on the uncompensated friction force. The higher is the friction, the higher acting force deviation |∆Fxmax | must be set up. We also compared the performance of the open-source solver lp solve and optimal Matlab 2015a solvers, namely linprog and quadprog, for LP (15), (16) and the equivalent QP, respectively, see Table 1. The implementation of ’simplex’ or ’dual-simplex’ algorithm with linprog does not provide an appropriate solution at all.
The optimization task is separated into two independent optimization programs, for LTV- and LTI-subsystems of (11), respectively. To minimize the absolute values of the predicted errors, the new variables 1 , 2 , 3 and ε1 , ε1 are introduced. These variables are positive and restricted by the corresponding absolute values of the chosen maximal allowed errors |∆xmax |, |∆Mbmax |, |∆Fxmax | and c max max |∆yh |, |∆y˙ h |, respectively. The considered system has the finite run time horizon NT and it must be taking into account. At first, we consider the case, when NT > N with the prediction horizon N and formulate the optimization task in the form of two independent linear programs, namely minpˆ1 1 + 2 + 3 subject to p1 , k) − xc,d (k) | ≤ 1 k = 1, .., N |xc (ˆ |Mb (ˆ p1 , k) − Mb,d (k) | ≤ 2 k = 1, .., N |Fx (ˆ p1 , k) − Fx,d (k) | ≤ 3 k = 1, .., N (15) 0 ≤ 1 ≤ |∆xmax | c max 0 ≤ 2 ≤ |∆Mb | 0 ≤ 3 ≤ |∆Fxmax | and minpˆ2 ε1 + ε2 subject to p2 , k) − yh,d (k) | ≤ ε1 k = 1, .., N |yh (ˆ |y˙ h (ˆ p2 , k) − y˙ h,d (k) | ≤ ε2 k = 1, .., N (16) 0 ≤ ε1 ≤ |∆yhmax | max 0 ≤ ε2 ≤ |∆y˙ h | with the optimization vectors T ˆ1 ˆ1 . . . h pˆ1 = 1 2 3 h , 5 N −1 T pˆ2 = [ ε1 ε2 yˆh,3 . . . yˆh,N −1 ] and planned trajectories (xc,d , Mb,d , Fx,d ), (yh,d , y˙ h,d ). ˆ 1 and y ˆ 0 depend on the initial state of The states h 0 SMC and can be calculated offline. When NT = N the prediction horizon N must be replaced by the decreasing prediction horizon N − n with the new actual step n = 0...N − 1. In general the initial value for the optimization vectors can be set zero. To achieve higher calculation efficiency, the optimization vectors are initialized by the
Optimization Tools quadprog linprog lp solve
Algorithm active-set active-set revised simplex
Calc. Time (s) ∼0.12 ∼0.06 <0.02
Table 1: Comparison of calculation time of optimal solvers implemented in Matlab 2015a CONCLUSION In this contribution, the tracking problem for the flexible SMC was solved by means of MPC. It was shown, that the simplified linear model of SMC can be used, to design an efficient model predictive control in the form of a linear program applying flatness-based parameterization. The parameterization of the system by the special flat output 5
Proceedings of the 9th MATHMOD 36 Vienna, Austria, February 21-23, 2018
XY - payload position
0.2 0
0.2
0.4 0.6 X (m) Mb - bending moment
Mb (Nm)
1
Trajectory Feedforward
0
0.5
1
1.5
2
20 0 -20
0.2
0.4
0.6 X (m) Mb - bending moment
0.8
1
Trajectory MPC
0
0.5
60 Trajectory Feedforward
40
0
0
-1
2.5
t (s) Fx - driving force
60
0.2
1
0
-1
Trajectory MPC
0.4
0
0.8
Mb (Nm)
0
Y (m)
Trajectory Feedforward
0.4
XY - payload position
0.6
Fx (N)
Y (m)
0.6
Fx (N)
A. Galkina et al. / IFAC PapersOnLine 51-2 (2018) 31–36
1
1.5 t (s) Fx - driving force
2
2.5
Trajectory MPC
40 20 0
0
0.5
1
1.5
2
2.5
t (s)
-20
0
0.5
1
1.5
2
2.5
t (s)
Fig. 3. Simulation results: left - approximate nonlinear model with feedforward control; right - approximate nonlinear model with MPC facilitates the formulation of an optimization task and provides the band structure of the corresponding Jacobian matrix. Such optimization task can be simply solved by means of any open-source linear optimal solver. The simulation of the designed MPC together with the nonlinear approximated model of SMC subjected to model uncertainties and external disturbances, such as uncompensated friction force, demonstrated a robust performance. The prediction horizon of MPC can be chosen high enough to avoid stability problems. For the future work, the implementation of MPC with the laboratory model of SMC and is planned. Also, the design of the optimal trajectory using the presented flatness-based parameterization will be investigated in the future.
Galkina, A., Gafur, I., and Schlacher, K. (2017). Model predictive control with linear programming for the strip infeed in hot rolling mills. In Proceedings of the 20th World Congress of the International Federation of Automatic Control, 1–6. Hajdu, S. and Gaspar, P. (2014). From modeling to robust control design of single-mast stacker cranes. Acta Polytechnica Hungarica, 11(10), 135–149. Luenberger, D. and Ye, Y. (2008). Linear and Nonlinear Programming. Springer, 233 Spring Street, New York, NY 10013, USA, 4 edition. Mayne, D., Rawlings, J., Rao, C., and Scokaert, P. (2000). Constrained model predictive control: Stability and optimality. Automatica, 36, 789–814. Rams, H., Sch¨oberl, M., and Schlacher, K. (2017). Optimal control of a single mast stacker crane. IEEE Transactions on Control System Technology. Shiobara, H. and Hori, N. (2008). Numerical exact discrete-time-model of linear time-varying systems. In Proceedings of International Conference on Control, Automation and Systems, 2314–2318. Staudecker, M., Schlacher, K., and Hansl, R. (2008). Passivity based control and time optimal trajectory planning of a single mast stacker crane. In Proceedings of the 17th IFAC World Congress. Zadeh, L. and Desoer, C. (1963). Linear System Theory: The State Space Approach. McGraw-Hill, New York.
REFERENCES Agrawal, A., Verschueren, R., Diamond, S., and Boyd, S. (2017). A rewriting system for convex optimization problems. https://arxiv.org/abs/1709.04494. Aschemann, H. and Schindele, D. (2010). Model predictive trajectory control for high-speed rack feeders. In In Model Predictive Control, 183–198. InTech. Bachmayer, M., Ulbrich, H., and Rudolph, J. (2011). Flatness-based control of a horizontally moving erected beam with a point mass. Mathematical and Computer Modelling of Dynamical Systems, 17(1), 49–69. Galkina, A., Gafur, I., and Schlacher, K. (2016). Modellpr¨ adiktive Regelung f¨ ur das Bandeinf¨ adeln in der Warmwalzstrasse. BHM Berg- und Huettenmaennische Monatshefte, 161(11), 520—-525. doi:10.1007/s00501016-0539-6. 6