Chance-constrained model predictive control for spacecraft rendezvous with disturbance estimation
- 格式:pdf
- 大小:1.23 MB
- 文档页数:12
Cross-coupled Model Predictive Control for Multi-axis Coordinated MotionSystemsXu Zhang1, 2, Kuanyi Zhu2, Xueyou Yang11. College of Precision Instrument and Optoelectronics Engineering, Tianjin University,Tianjin 300072, ChinaE-mail: redheartzx@2. School of Electrical and Electronic Engineering, Nanyang Technological University,Singapore 639798, SingaporeAbstractA cross-coupled model predictive controller for multi-axis coordinated motion systems is presented by incorporating cross-coupled strategy into the model predictive control architecture. The control law is generated by minimizing a cost function for the augmented system model in which the synchronization errors are embedded. Simulation and experimental results conducted on a high-precision positioning system with two permanent magnet linear motors (PMLMs) demonstrate the effectiveness of the proposed approach.1.IntroductionHigh-precision multi-axis motion systems are essential elements in modern advanced manufacturing. The term “axis of motion” refers to one degree of freedom, or forward and backward motion along one direction, which may be linear or rotary motion. When two or more axes of motion are involved on a single machine, that machine is employing multi-axis motion. In precision motion control, besides good tracking performance of each individual axis is absolutely necessary, synchronization or coordination between multiple axes becomes another challenging problem needed to be solved [1].Demand for high-precision control of multi-axis motion systems calls for the development of high performance positioning devices and accompanying control approaches. With the attractive features such as high dynamic performance and most importantly high positioning precision, permanent magnet linear motors (PM LM s) are increasingly widely being applied in multi-axis motion systems [2]. To harness the potential of PM LM s, the main challenge lies in the control approach.Cross-coupled control is an appropriate strategy used for coordinated control of processes with multiple motion axes. This concept has been successfully applied to reducing contour error of multiple motion axes in machine tool control [3]–[6], minimizing differential velocity error of two driving wheels in a mobile vehicle [7], and synchronizing the multi-robot for assembly tasks respectively [8]. However, the principal shortcoming of the existing control techniques is that they cannot explicitly incorporate plant model uncertainties in the formulation of the high precision control system when using the traditional PID or fussy method.M odel predictive control (M PC) algorithm is an effective approach to improve motion tracking and coordination performance of multi-axis systems due to its control acumen to deal with the multivariable constrained control problems. In the past two decades, the theory of MPC has been well developed nearly in all aspects, such as stability, nonlinearity, and robustness [9]–[13]. Several methods based on M PC are presented for better positioning and tracking precision of a single axis or motor [14]–[17]. A cross-coupling framework of generalized predictive control is investigated by compensating both tracking error and synchronous error between the process and the reference model [18]. In this paper, design and implementation of a cross-coupled model predictive control algorithm are presented for multi-axis motion systems. This paper is organized as follows. In Section 2, the dynamics of the PM LM is presented. Process model and MPC is formulated in Section 3. In Section 4, synchronization strategy is proposed. The cross coupled M PC controller for multi-axis coordinated motion systems is presented in Section 5. Simulations and experiments are performed in Section 6, and finally, conclusions are given in Section 7.International Conference on Advanced Computer Control2. Dynamics of the PMLMThe dynamics of the PM LM can be viewed as comprising of a dominantly linear model and a nonlinear remnant [19]. In the dominant linear model, the mechanical and electrical dynamics of a PM LM can be expressed as follows:d m MxDx F F ++= (1) ()ae a a a dI K xL R I u t dt++= (2)m f a F K I = (3)where x(t) denotes motor position; M ,D ,F m , and F d denote mechanical parameters, i.e., mass of the moving thrust block, viscosity constant, generated force, and system disturbance, respectively; u(t),I a ,R a , and L a denote electrical parameters, i.e., time-varying motor terminal voltage, armature current, armature resistance, and armature inductance, respectively; K f denotes an electrical–mechanical energy conversion constant; and K e is back electromotive force (EM F) voltage. The delay due to electrical transient response may be ignored since the electrical time constant is typically much smaller than the mechanical one. The following is the simplified model:()()()x axt bu t t ξ=++ (4) wheref e a K K R D a RM+=−fa Kb R M =()dF t M ξ=−.Note that Ӷ(t) includes extraneous nonlinear effectsthat may be present in the physical structure. Among them, the prominent nonlinear effect is due to ripple force arising from the magnetic structure of the PMLM and other physical imperfections, which can be viewed as a sinusoidal function of motor position with a period of Ԁ and amplitudes of A r and B r , i.e.,cos()sin()ripple r r F A x B x ωω=+ (5)where F ripple is ripple force. In particular, ripple force ismore complex in shape, e.g., due to variations in the magnet dimension.3. Process model and MPCA process with multiple PM LM s can be describedby the linear discrete-time CARIMA model:111()()()(1)()()/i i i i i i A q x t B q u t C q t ξ−−−=−+Δ(6)where t is the normalized discrete time, i denotes the corresponding motor of the process, x i (t) and u i (t) are the position and control input of the motor, respectively, Ӷi (k) represents system disturbance andunmodeled dynamics, A i (q -1)and B i (q -1) are the systempolynomials in the backward shift operator q -1,C i (q -1) is a polynomial used for robustness of the controller, and Ӕ is the differencing operator.In model predictive control, minimization of a cost function yields the control law. The cost function is:211{(()(|))Hpn c i i i k I w t k x t k t ===+−+¦¦,21((1|)}Hci i i k u t k t λ=+Δ+−¦ (7)where w i (t + k) is the reference signal, H p and H c are the prediction and control horizons, respectively, ӳi is a non-negative weighing factor for the control input. Based on MPC technique, the future output x i (t + k|t)can be expressed as(|)(1|)()i i ik ik x t k t G u t k t p t +=+−+,1,2,3,......,p k N = (8)where11k ik ik ik G G g q −+−=+ (9)()()(1)i ik ik ik i iiF H p t x t u t C C =+Δ− (10)G ik , F ik and H ik satisfy the Diophantine equationsk i i ik ik C A E q F −=Δ+ (11) k i ik i ik ik B E C G q H −=+ (12)with the corresponding degrees() () 1() max((), ())() max((), ())1ik ik iki i ik i i E G k F A C k H B C δδδδδδδδ==−=−=−°®°¯(13) The matrix version of the predictor (8) is obtained asi i i i X G U P =+ (14)where[( 1|) ( 2|) ... ( |)]T i i i i p X x t t x t t x t N t =+++ (15) [( |) ( 1|) ... ( 1|)]T i i i i u U u t t u t t u t N t =ΔΔ+Δ+− (16)12 [() () ... ()]p N T i i i i P p t p t p t = (17)0100(1)(2)(1)(2)()000u u p p p u i i i i N i N i i N i N i N N i g g g G g g g g g g −−−−−…=…ªº«»«»«»«»«»«»«»«»¬¼%##%%%#%%# (18)4. Synchronization strategyIn multi-axis coordinated motion systems, both tracking performance and synchronization performance are specified. However, if the motors follow different trajectories, i.e., r i (t) r j (t), then a scalar of r i (t)/r j (t) should be introduced to describe the relationship between the ith motor and the jth motor accurately. So the synchronization error can be defined as() (() ()) ij ij i ij j t x t x t φαβ=− (19) where өij and Ӫij are the coupling and synchronizationfactors between the ith motor and the jth motor, respectively. M oreover, even if Ĩij (t) = 0 may be achieved at an instant time, the overall system may not be synchronized perfectly due to ӔĨij (t + 1) 0,where ӔĨij (t + 1) = Ĩij (t + 1) íĨij (t). In this case, the difference between the consecutive synchronization errors should also be considered as1() (1 )()ij ij t q t ψφ−=− (20)denoting[( 1|) ( 2|) ... ( |)]T ij ij ij ij p t t t t t N t φφφΦ=+++ (21) [( 1|) ( 2|) ... ( |)]T ij ij ij ij p t t t t t N t ψψψΨ=+++ (22)For the synchronization error and its difference, bysetting ӿij(t + 1|t) = 0,it can be obtained thatij ij L Ψ=Φ (23)where00001100011100011L −=−−−ªº«»«»«»«»«»«»¬¼""%##%%" (24)The cost function is introduced as1{()()}i i i i i i i nT T c i W X W X U U I λ=−−+=¦111n nij ij i j i T −==++ΦΦ¦¦ȁ (25)where[( 1) ( 2) ... ( )]T i i i i p W w t w t w t N =+++ (26)TI L L υη=+ȁ (27)and I is the unit matrix.5. Cross-coupled MPC controllerBy minimizing the cost function (25), the Cross-coupled MPC control law can be obtained:()()21111122122222312n m m m m n m m m m n n n W P P P U U W P P P U W P αβαβ=−==ΩΓ§ªº−−Λ−ªº¨«»«»¨«»−−Λ−«»¨«»«»¨«»«»¨«»¨¬¼−¬¼©¦¦##()()12222211210m m m m m n mn mn mmn n m P P P P αββαββ=−=·ªº¸«»Λ−¸«»+¸«»¸«»¸«»¸Λ−¬¼¹¦¦# (28)where Ө,ӓ,Өij are expressed as follows1111n n nn ΩΩΩ=ΩΩªº«»«»«»¬¼"#%#" (29)100T Tn G G Γ=ªº«»«»«»¬¼"#%#" (30)2, ,Ti i ij Tij ij j i i G G I if i jG G if i j λαρΟ+=ΩΛ<°=®−°¯ (31) 122211()i nmi mi im m m i I αβα−==+Ο=+Λ+¦¦ (32)The block diagram shown in Figure 1 illustrates thestructure of cross-coupled MPC controller.Figure 1. Structure of cross-coupled MPCcontroller6. Simulation and experimentFigure 2. High-precision positioning controlsystemThe effectiveness of the proposed control approach will be illustrated by computer simulations and experiments on a high-precision positioning system with two PMLMs, as shown in Figure 2. The dynamics of the positioning system can be described by11120.653()11.006()x xt u t =−+ ()11.418sin 20.074x π+ ()10.313cos 20.074x π− (33)2225.917() 5.048()x xt u t =−+ ()20.676sin 20.051x π+ ()2cos 0.19120.051x π− (34)With a sampling interval T s = 0.0008 seconds, the corresponding linear discrete-time system models are given by111()() ()(1) ()/(1),i i i i i i A q x t B q u t C t q ξ−−−=−+− 1, 2i = (35)where11211122111112112() 1.000 1.9840.984() 1.000 1.9950.995() 3.502 3.483() 1.613 1.610 10.8,,,,.A q q q A q q qB q q B q qC C q −−−−−−−−−−−=−+=−+=−=−==−In this study, two PM LM s are desired to track thesame square wave. The parameters of the cross-coupled M PC controller are set as: N p = 10; N c = 4; ӳ=50,ӵ = 20,ӯ = 10,ө12 = 1,Ӫ12 = 1. Figure 3-5 show the desired tracking trajectories of two PMLMs, the simulation and experimental results for the cross-coupled MPC controller, respectively.Figure 3. Desired tracking trajectoriesFigure 4. Simulation results for the cross-coupled MPC controllerFigure 5. Experimental results for the cross-coupled MPC controller7. ConclusionTo achieve better synchronization performance while maintaining good tracking response, a cross-coupled M PC controller for multi-axis coordinated motion systems has been presented. Computer simulations and experiments on a high-precision positioning system with two PMLMs have shown that the improved performance of the control system can be obtained with the control parameters to be set as follows: first, N p ,N c and ӳare selected appropriately so that the control system provides a stable andsatisfactory tracking response and then, the weighting factorsӵandӯare introduced to reduce the synchronization errors caused by the mismatched dynamics, system disturbances, and so on. Although the same desired trajectory is considered in this study, different trajectories can be realized easily by setting the synchronization factorsөand Ӫʳappropriately.8. References[1] T.C. Chiu and M. Tomizuka, “Coordinated position control of multi-axis mechanical system”, Trans. ASME, J. Dyn. Syst. Meas. Control, 1998, vol. 120, no.3, pp.389–393. [2] K.K. Tan, T.H. Lee, H.F. Dou, and S.N. Huang, Precision Motion Control: Design and Implementation (Advances in Industrial Control). London, U.K.: Springer-Verlag, 2001.[3] J. Borenstein and Y. Koren, “Motion control analysis ofa mobile robot”, Trans. ASME, J. Dyn. Syst. Meas. Control, Jun. 1987, vol. 109, no. 2, pp. 73–79.[4] Y. Koren and C.C. Lo, “Variable-gain cross-coupling controller for contouring,” Ann. CIRP, 1991, vol. 40, no. 1, pp. 371–374.[5] S.S. Yeh and P.L. Hsu, “Analysis and design of integrated control for multi-axis motion systems” , IEEE Trans. Control Syst. Technol., May 2003, vol. 11, no. 3, pp. 375–382.[6] P.K. Kulkarni and K. Srinivasan, “Cross-coupled control of biaxial feed drive servomechanisms”, Trans. ASME, J. Dyn. Syst. Meas. Control, 1990, vol. 112, no. 2, pp. 225–232.[7] L. Feng, Y. Koren, and J. Borenstein, “Cross-coupling motion controller for mobile robots,” IEEE Control Syst. Mag., Dec. 1993, vol. 13, no. 6, pp. 35–43.[8] D. Sun and J.K. M ills, “Adaptive synchronized control for coordination of multirobot assembly tasks,” IEEE Trans. Robot. Autom., Aug. 2002, vol. 18, no. 4, pp. 498–510. [9] D. W. Clarke, C. Mohtadi, and P. S. Tuffs, “Generalized predictive control—Part I. The basic algorithm,” Automatica, Mar. 1987, vol. 23, no. 2, pp. 137–148.[10]D. W. Clarke, C. Mohtadi, and P. S. Tuffs, “Generalized predictive control—Part II. Extensions and interpretations,” Automatica, Mar. 1987, vol. 23, no. 2, pp. 149–160.[11] M. Morari and J. H. Lee, “Model predictive control: Past, present and future,” Comput. Chem. Eng., 1999, vol. 23, no. 4/5, pp. 667–682.[12] D. Q. Mayne, J. B. Rawlings, C. V. Rao, and P. O. M. Scokaert, “Constrained model predictive control: Stability and optimality,” Automatica, Jun. 2000, vol. 36, no. 6, pp. 789–814.[13] J. M. M aciejowski, Predictive Control: With Constraints. London, U.K.: Prentice-Hall, 2002.[14] L. Zhang, H.R. Norman, and W. Shepherd, “Long range predictive control of a current regulated PWM for induction drives using the synchronous reference frame,” IEEE Trans. Control Syst. Technol., Jan. 1997, vol. 5, pp. 119–126.[15] K.S. Low, K.Y. Chiun, and K.V. Ling, “Evaluating Generalized Predictive Control for a Brushless dc Drive,” IEEE Trans. Power Electron., Nov. 1998, vol. 13, pp. 1191–1198.[16] K. Hentabli, M.E.H. Benbouzid, and D. Pinchon, “M ultivariable state-space CGPC application to induction motor control,” in Proc IECON’97, 1997, vol. 1, pp. 181–186.[17] M. Kinnaert, “Adaptive generalized predictive controller for MIMO systems,” Int. J. Control, 1989, vol. 50, no. 1, pp. 161–172.[18] K.Y. Zhu and B.P. Chen, “Cross-coupling design of generalized predictive control with reference models,” in Proc. Inst. Mech. Eng., 2001, vol. 215, pt. I, no. 4, pp. 375–384.。
具有输出限制的纯反馈系统的神经网络控制史昱;尹丽子【摘要】为了研究一类具有输出限制的不确定非线性纯反馈系统的自适应神经网络追踪控制问题,利用神经网络的非线性逼近能力与自适应控制的反推法给出该系统的自适应控制器;利用障碍Lyapunov函数与隐函数存在定理进行控制器的设计.结果表明,该控制方法保证了闭环系统所有信号的半全局一致最终有界性.%To investigate the adaptive neural network tracting control problem of a class of uncertain nonlinear pure-feedback systems with output constraints, an adaptive controller of the systems was provided by using the ability of Neural Network approximation and the adaptive backstepping techniques.The controller was designed by the barrier Lyapunov function and the implicit function theorem.The results show that the developed control scheme guarantees semiglobally uniform ultimate boundedness of all the signals in the closed-loop systems.【期刊名称】《济南大学学报(自然科学版)》【年(卷),期】2017(031)005【总页数】7页(P394-400)【关键词】纯反馈系统;输出限制;障碍Lyapunov函数;隐函数定理;自适应控制【作者】史昱;尹丽子【作者单位】山东交通学院理学院,山东济南250357;济南大学数学科学学院,山东济南250022【正文语种】中文【中图分类】O231.2多层神经网络、径向基(RBF)神经网络、高阶神经网络可以在紧集内以任意精度逼近非线性函数,具有良好的函数逼近能力,被广泛应用于系统函数未知的控制工程[1-2]。
具有可参数化不确定性系统的对偶自适应模型预测控制曹文祺; 李少远【期刊名称】《《控制理论与应用》》【年(卷),期】2019(036)008【总页数】10页(P1197-1206)【关键词】模型预测控制; 自适应控制; 对偶控制; 不确定性; 动态系统【作者】曹文祺; 李少远【作者单位】上海交通大学自动化系; 系统控制与信息处理教育部重点实验室上海200240【正文语种】中文1 引言模型预测控制(model predictive control,MPC)可以有效解决约束系统的控制问题,所以被广泛研究、运用[1].在这个方法中,人们使用对象模型预测系统未来的表现,从而得到动态系统的优化控制.因此模型质量对MPC应用的效果有直接的影响,建模过程耗时又复杂.然而,利用模型准确描述被控对象的动态过程并不简单:往往难以获取其实际结构、准确参数,而真实对象中还会存在难以描述的高频成分,即不可建模动态;上述因素使基于模型的MPC控制系统具有未来不确定性.尽管MPC自身由于其滚动时域特性,能保证对充分小不确定性在一定程度上的鲁棒性,但这并不够满足实际中的需求.在MPC领域,对闭环系统的不确定性处理,是理论研究与应用技术面临的共性问题[2].针对带有不确定性的过程对象,鲁棒模型预测控制(robust model predictive control,RMPC)一直被广泛研究.鲁棒控制的思想是:在算法设计初期就将不确定性的最大界考虑进去,得到能使控制系统在一定范围内的扰动下依旧稳定运行的控制律.已经存在太多有关RMPC的研究,可参见文献[3-5],本文不展开具体讨论.然而鲁棒控制的性能,总会受到固定不变的初始模型的质量限制[6],具有较强保守性,在模型发生较大变化时当前的控制方法就不再适用.利用自适应工具对系统在线更新是解决不确定性问题的另一种方法.迄今,自适应MPC的研究仅受到很少关注[7].自适应控制利用辨识技术,主动修正控制器本身的参数和策略,在原始模型未知或不准、控制环境改变时仍能得到令人满意的控制效果.传统自适应MPC方法将辨识和控制分离,没有在控制中考虑系统的不确定性,特别是未建模动态带来的影响.Feldbaum[8]于1960年提出,对系统的信号控制是一个具有对偶特性的过程:一方面尽可能地使控制结果更接近目标要求;另一方面利用在控制过程中捕获(probe)的信息优化控制器,从而提升控制器的未来表现.这个信息就是系统的不确定性.实际上,这两个方面互相冲突:前者以减小系统不确定性为目标,后者要增加对不确定性的利用;兼顾两者的控制即为“对偶控制(dual control)”.近10年对偶控制才被引入MPC的研究中.早期的工作比如Marafioti等人[9]以及ˇZ´aˇcekovˇa等人[10]通过持续激励输入信号实现MPC的对偶特性.Kim和Sugie[6]在前人RMPC的基础上,对参数估计误差建立了一个单调递减的边界数值函数,并给出了一个减小参数不确定性的、非传统形式的自适应的鲁棒控制方法.王超等人[11]基于状态反馈Tube不变集RMPC算法,利用自适应集元滤波在线估计系统过程噪声边界及状态可达集,给出了确保系统状态鲁棒渐近稳定,并收敛于终端干扰不变集的鲁棒自适应MPC算法.最近的工作中,人们通过修改目标函数来利用系统的不确定性.Houska等人[12]针对含状态估计器的MPC控制,提出一种名为self-reflective的非线性MPC控制策略,通过在目标函数中加入未来状态估计方差的函数实现对偶控制;Feng等人[13]延续了这个工作,为self-reflective MPC提供了目标函数权重自动确定的实时方法.前人有关对偶MPC的研究,基于固定的对象模型,虽然利用了系统中存在的不确定性,但当对象模型不准确或发生改变时,不确定性会增大而影响控制.持续激励输入的方法,无法充分描述系统当前的不确定性,有时还会因过多信息的引入为控制带来更多不便.Heirung等人[14]于2015年率先研究了对偶自适应MPC(dual adaptive MPC)控制:针对参数估计不确定性,基于确定性等价(certainty-equivalence,CE)[15]自适应MPC,在其目标函数中加入未来参数估计的信息矩阵或协方差矩阵的函数项,优化性能指标的同时减小自适应MPC中的参数不确定性.这个开创性的工作利用自适应控制,解决了对象参数不确定性对控制性能带来的影响.在控制中利用不确定信息的对偶控制环节,虽然提高了未来辨识的准确性,但并没有对自适应控制的缺点进行补偿,即无法解决一直存在的未建模动态不确定性问题.在2017年的发表的另一论文中,Heirung等人[16]针对随机系统中存在的不确定性问题,基于正交基模型,通过添加对未来跟踪误差的统计描述项,建立了一个以最小化期望性能损失为目标的随机最优控制问题,并将其近似为二次约束二次规划(quadratically constrained quadratic program,QCQP)问题求解.将不确定性相关项加入目标函数,诚然实现了不确定性的减小;但是原目标函数中的性能指标和不确定性相关项存在耦合,即使可以调整权重,目标函数的最小化依旧不能明确体现性能指标或不确定性自身单独的优化程度.在实际应用中,修改目标函数的做法也可能无法满足人们对性能指标较为严格的要求.本文提出一种减小系统中可参数化不确定性影响的对偶自适应MPC算法:使被控对象的数学模型尽可能准确的同时,当前MPC控制系统的未来不确定性也在减小.具体的做法是:首先将动态过程中的不确定性参数化表达;然后描述在自适应MPC 的迭代过程中对不确定性的边界要求,得到约束不确定性的关于时间(迭代步数)的边界函数;最后在每一步MPC控制中,将当前时刻对不确定性的边界要求作为额外约束添加其中,得到在参数估计准确,且系统不确定性在用户要求范围内的最优控制输入.该算法结合了鲁棒控制与自适应控制的策略,兼具两者的优点,也弥补各自的缺点,均衡考虑了模型质量和未建模动态的影响.不确定性被参数化表达;不同于他人利用不确定信息的方法,添加不等式约束的形式,在控制中实现控制性能优化与不确定性信息优化的分离,从而得到对目标性能的更为准确的控制.值得一提的是:约束参数化不确定性的边界函数,可由用户自由设计,以满足对减小不确定性影响的不同的具体需求;几种边界函数的设计建议也将被给出.本文进行了MATLAB仿真实验与敏感性分析,在与传统的确定性等价自适应MPC的对比实验中,对偶MPC算法在控制过程中会得到更为平滑的输入、输出,且具有更好收敛性.本文结构如下:第2章讨论对偶控制的原理,介绍被控对象的数学模型及确定性等价自适应MPC 的实现流程;第3章分别就自适应模型预测控制中不确定性的参数化方法、算法的具体实现以及边界函数的设计介绍本文的自适应对偶MPC 算法;第4 章是MATLAB仿真对比实验和敏感性仿真实验的结果及分析;第5章总结并展望未来可进行的研究工作.附录中是定理1的证明.2 背景及问题描述2.1 对偶控制的原理本文面向具有未建模动态及被控对象未知参数的控制系统,故需对被控对象进行参数估计.图1为传统的确定性等价自适应MPC的结构框图.其中,参数估计模块根据历史控制输入和系统输出,用辨识方法计算参数的最优估计值;CE MPC控制模块认为参数估计模块估计的参数值就是被控对象模型的真实参数值,对当前所谓的“参数真实的对象”进行MPC控制.将控制后得到的新的输入输出数据传入参数估计模块,更新模型中的参数估计值,再进行MPC控制.如此循环往复,即为确定性等价自适应MPC的基本流程.显然,图中由虚线标注的部分存在未来参数、未建模动态等不确定性相关信息,但是传统的控制方法并没有对不确定性相关信息加以利用.图1 确定性等价自适应MPC的结构框图Fig.1 Block diagram illustrating a certainty-equivalence MPC structure图2为对偶自适应MPC的结构框图.对比确定性等价自适应MPC结构框图,对偶MPC(dual MPC)控制模块利用了虚线框中存在的不确定性相关信息,得到可以使未来不确定性影响减小的控制输入.更具体地说,除了完成普通MPC的工作,dual MPC控制模块还利用当前系统中可以获取的数据,去描述未来系统中存在的不确定性,在计算控制律的同时减小参数化的不确定性大小,最终得到可以减小未来不确定性影响的控制律.图2 对偶自适应MPC的结构框图Fig.2 Block diagram illustrating a dual MPC structure2.2 被控对象的描述本文考虑的是离散的线性时不变单输入单输出(single input single output,SISO)系统.根据上述假设,将其模型描述如下:其中:t=1,2,3,···为离散时间点,y(t),u(t),v(t)分别为被控对象的输出、控制输入和t 时刻的不可测噪声(未建模动态);系数和(b1/=0)为未知的参数;假设v(t)为独立同分布的高斯随机变量,均值为0,方差为r,即下面将模型简写为分别包含历史输入输出数据和未知参数.2.3 确定性等价自适应MPC下面分别对图1所示的确定性等价自适应MPC的参数估计部分和CE MPC控制部分的算法进行介绍.本文选择的辨识方法为递推最小二乘法.2.3.1 递推最小二乘法使用最小二乘法估计参数可以实现统计意义上的最优.辨识本文被控对象的参数时,目标函数为该目标函数可以通过求导直接得到最小值.但是在线控制的过程中,输入输出数据量逐渐增加,为了减少计算的不便,选择利用新输入输出修正上一步辨识结果以得到新结果的递推方法.令Φ(t)=[φT(0)φT(1)··· φT(t)]T,记为t(t=1,2,···)时刻最小二乘法得到的对参数的最优估计,通过推导易知t时刻参数估计的协方差矩阵可表示为定义矩阵易得递推公式其中I为与P(t)同阶的单位矩阵.使用递推最小二乘法得到了模型参数后,系统往往会因为无法描述随机干扰而按如下模型被使用:2.3.2 确定性等价MPCMPC一般会有类似如下形式的二次有限时域目标函数:其中:N为预测时域(式(10)取控制时域等于预测时域),w1,w2为权重参数,φ(k)包含对i=k-na+1,k-na+2,···,k,在i ≤t 时的记录输出y(i),i >t时由模型计算得的预测输出y(i|t);以及对j=k-nb+1,k-nb+2,···,k,j <t时的记录输入u(j)和j ≥t时的自由变量u(j).需要说明的是,目标函数可根据需要进行调整,如:取不相等的预测时域与控制时域,为输出、输入变量设置非零参考值,添加控制输入变化量有关的项(有时也可通过将模型转化为增量模型,直接用u(t)表示控制输入变化量)等.在确定性等价MPC中,t时刻参数值为辨识方法估计所得的(t),因此式(10)中的θ应选用(t).记k时刻MPC控制中,由(t)计算得到控制输出为(k+1|t),含由(t)计算所得预测输出或输入变量的回归向量为(k).又,大部分MPC控制均包含对输入、输出的上下界约束;故确定性等价MPC的优化问题可被描述为其中umin,umax,ymin,ymax分别为输入与输出的上、下界.式(12)是一个凸二次规划(quadratic programming,QP)问题.问题共有线性约束(3na+2nb)N个.其中:等式约束个数为naN,不等式约束个数2(na+nb)N.解该问题,得到最优开环控制序列取序列的第一个值uCE(t),令u(t)=uCE(t),对未知系统进行控制,得到新的系统输出;再利用新获得的数据更新参数估计值,如此循环往复,构成完整的确定性等价MPC控制.3 对偶自适应MPC本文具有可参数化不确定性系统的对偶自适应MPC算法通过对由未建模动态等引起的系统未来不确定性影响进行数学描述,在自适应MPC的控制环节添加了有关不确定性的约束.这样做保证被控对象参数估计准确的同时,实现了具有鲁棒性能的MPC控制,弥补了传统自适应MPC未充分考虑未建模动态的缺陷;最终得到未来不确定性影响满足一定范围的使控制性能指标函数最优的控制律,从整体上提升了控制性能.由于主动对未建模动态不确定性进行约束,某种程度上增加了未来模型估计的准确性,所以本文算法也可使优化问题中由计算的预测值描述的性能指标更接近下一时刻性能指标的真实值,同时提高了系统未来真实输出满足式(12d)中对预测输出的上下界约束的可能.定义t时刻在MPC问题中添加的不确定性约束为其中:ρ(t+1)为t 时刻描述未来不确定性的函数,与未来不确定性成正相关;ρub(t+1)为根据用户需求预定义的边界函数,代表当前时刻用户能够允许的未来不确定性的最大值.第3.3节将对边界函数的设计进行讨论并给出几种可行的方法.3.1 不确定性的讨论为了利用控制中的不确定性,首先需要尝试将其对系统控制的影响参数化表达.本文欲在MPC 控制环节对系统未来的未建模动态添加约束,实际上就是在t时刻想办法对式(3)中的随机噪声v(i)(i=t,t+1,···)进行具体描述.在真实系统中,当参数准确时,式(3)反映了回归向量(历史输入、输出数据)和输出数据之间的函数关系;如果已知回归向量和输出数据,将数据看作参数,则式(3)反映了随机干扰v(t)与原模型参数θ的计算关系:对一对已知的φ(t-1),y(t)数据,若已知θ的大小,则可求此时v(t)大小.在辨识过程中用到的模型其实是式(9),经过参数辨识得到的最优(t)值,使回归向量和输出数据尽可能满足由(t)描述的式(9)中的模型.观察式(3)(9)的差别,可以发现(t)中其实包含了对随机干扰的描述;所以不难理解,参数真值θ与最优估计值(t)的差别从一定程度反映了随机干扰,也就是未建模动态的大小;未来时刻的参数估计值(t)与真值θ的偏离程度,体现了未来系统中不确定性的大小.自适应控制中,参数估计的方差矩阵P(t)描述估计值(t)与真值θ 的偏离程度.因此可以说,未来时刻的参数估计方差矩阵正相关地描述了未来系统中的不确定性大小:P(t)越大,(t)偏离θ真值越大,系统受到未建模动态的影响就越大,未建模动态不确定性就越大.对t时刻的MPC问题,取未来Ne(Ne∈N∗)个时刻的参数估计方差矩阵的函数组合构成ρ(t+1),来描述未建模动态引起的不确定性:其中:f(P(k))为P(k)的函数(需满足f(P(k))能正相关地描述参数估计的不确定性或P(k)的大小),δ ∈(0,1]为f(P(k))的遗忘因子.使用遗忘因子,可以弱化遥远未来影响,增强被描述得更加准确的近期未来的影响.本文建议Ne的取值小于预测时域N,过大的Ne会因模型不准、控制的主要目的是满足控制性能指标等原因而没有意义. 将式(14)代入式(13)便能得到一个完整的有关未来不确定性的约束.第3.2节中会给出f(P(k))的设计建议并在后文及仿真中选用其中一种设计方法.3.2 对偶MPC优化问题下面给出本文对偶自适应MPC算法的MPC优化问题.为了将有关未来参数估计协方差矩阵的不等式约束引入优化问题,需向优化问题中加入有关P(k)(k=t+1,···,t+Ne)计算的等式约束.由式(12)-(14),以及方差矩阵的迭代公式(7)(8a),具有对偶特性的MPC优化问题如下:约束于其中需要注意的是:P(k+1)(k=t,t+1,···,t+Ne-1)不是k+1时刻协方差矩阵的真实值,而是由自由变量u(j)(j ≥t)、历史输入输出数据,以及由式(15b)计算的预测输出共同计算表示的变量;可以较为接近地表示未来协方差矩阵的值.这个对偶MPC问题是一个非线性规划(nonlinear programming,NLP)问题;对比确定性等价MPC的QP优化问题,本文算法增加了对方差矩阵的利用,为优化问题添加了3组非线性约束,很大程度增加了问题的复杂程度.定义辅助变量z(k)∈R(na+nb)×1和ζ(k)∈R1 如下:又考虑到方差矩阵P(k+1)为对称矩阵的性质,K(k+1),P(k+1)的迭代公式可化为简单的二次约束形式:f(P(k))有多种可行的设计方法,理论上,对本文系统,能正相关地描述不确定性或P(k)大小的设计均可行.行列式、各种范数、矩阵的迹等都可以用来衡量矩阵的大小[14],显著的不同之处在于计算的复杂程度.矩阵的范数涉及代价高昂的非线性计算,增加了原本就因P(k)的存在而颇受影响的计算复杂程度.行列式提供了矩阵范数的一种替代.但是当阶数超过三时,行列式不再能由简单的公式直接计算;常规计算方法需要展开行列式,求元素与代数余子式乘积之和,包含大量乘法,十分复杂.一个简单的行列式计算方法是求所有特征值的乘积,但是特征值的计算也是非线性且复杂的.可以通过利用矩阵指数公式使用行列式:其中det(M),exp(M),tr(M)分别表示矩阵M的行列式、矩阵的自然指数和矩阵的迹.矩阵的迹可以直接通过对矩阵对角线所有元素求和得到,易于计算,且该计算过程无需向优化问题中添加额外约束.采用式(18)是一种可行且计算较为简便的f(P(k+1))的选择.本文直接选则矩阵的迹作为f函数:如此便可通过计算P(k+1)对角线元素之和得到函数值.向优化问题(15)中代入式(16)-(17)(19),得到重写的优化问题:至此已完成对偶MPC优化问题的全部推导.对比确定性等价MPC的优化问题(12),对偶MPC的优化问题中增加了二次等式约束共(na+nb+1)2Ne项;因为矩阵P(k+1)对称,故可以只取其左下三角或右上三角(含对角线)部分的元素进行约束,二次等式约束减少为[(na+nb)2+5(na+nb)+2]Ne/2项.在后添加的未来不确定性的不等式约束项(20i)中,ρub(t+1)已定义,可认为是常数;不等式左端可认为是不同时刻未来方差矩阵的对角线元素的线性组合;故式(20i)为线性不等式约束.经过辅助变量的引入和函数f(P(k+1))的选取,本文对偶MPC的初始NLP问题(15)被化为一个二次约束二次规划(quadratically-constrained quadratic pro gramming,QCQP)问题.现在已有求此类问题全局最优解的算法,如GloMIQO算法[17].3.3 自定义边界函数的设计本节给出边界函数ρub(t+1)的设计建议.因为ρub(t+1)约束系统参数化的不确定性,所以一种直观且简单的选取就是取常数ρc(ρc>0),则不确定性约束(20i)变为式(21)可以实现在控制的过程中,系统参数化的不确定性一直小于ρc.还可以将边界函数设置为分段函数:其中:ρ1>ρ2>0为设定的两个常数,k=1为自适应MPC的初始时刻,tm为设定的转换时刻;因为辨识的初期辨识参数不够准确,可能存在较大的偏差,所以这样设置能够在自适应MPC控制的前期中更加充分地利用系统信息,实现更快的辨识;随着辨识接近收敛,再将不确定性的边界设置得更小,以实现未来控制的良好的鲁棒性. 分析本文中的不确定性约束(20i),定义函数则有如下定理:定理1 当P(k+1)(k=t,t+1,···,t+Ne-1)为由式(15e)(15f)描述的、k+1时刻对辨识方差矩阵的估计时,式(23)中定义的函数g(t)为单调递减函数.定理1的证明见附录.由定理1可知,当不确定性参数化为式(14)(19)所示的未来时刻方差矩阵的迹的加权和时,式(21)约束中的常数ρc其实只对初始时刻的g(t)进行了约束,因为当k >1 时,必有g(k)<g(1),又通过控制,有g(k)<g(1)<ρc.所以不妨将边界函数也定义成一个单调递减的函数,从而实现对每一次迭代的不确定性都进行一定程度的约束,最简单的定义方式比如使用比例因子λ ∈(0,1],取也可在一段时间后停止对不确定性的约束,或当不确定性函数值达到一定指标时停止对其的进一步控制,如令除上述列出的边界函数外,还可设计很多其他更为精准的、可满足用户不同需求的函数,本文不再展开讨论.4 仿真及结果4.1 可容许性及可行性分析参数估计中的可容许性(admissibility)问题运用在自适应控制时[18],表明参数空间中存在一些辨识出来的模型参数,会导致不良的控制性能,其中一个典型的例子就是零极点相消.本文明显不可容许的一组参数估计就是=.对式(1)模型,当y(t-i)(i=1,2,···,na)的系数均为0时,系统将不再表现出动态特性;若u(t-j)(j=1,2,···,nb)的系数均为0,系统变为不可控系统;此时MPC控制器将不会给出非零的控制输入.所以初始化参数估计时应避免取零向量.假设用户对控制输入、控制输出设定的上下界约束可以使确定性等价自适应MPC 问题(12)具有可行域.分析本文添加了参数化不确定性边界约束的对偶自适应MPC控制,在自适应控制的过程中,不确定性函数会随着辨识的收敛而越来越小,对因此造成的控制输入的约束也会在控制过程中逐渐放宽;可以认为:只要在自适应控制的初始,不选择过分严格的边界约束值,则本文的对偶MPC控制问题可行.4.2 仿真实验本文对如下具体系统进行仿真实验:其中a/=0且b/=0.则有向量t时刻对t+1时刻输出的预测为除非特殊声明,所有仿真取a=1.1(不稳定系统),b=1,噪声的方差r=0.01.取y(0)=1,参数估计的初值为(0)=-0.1,(0)=-0.1,递推最小二乘法中取初始值P(0)=10-3I.在MPC 目标函数中取ω1=ω2=1,N=10,令输入、输出的上下界分别为umin=-5,umax=5,ymin=-10,ymax=10.首先进行确定性等价自适应MPC与本文算法的对比仿真实验.对对偶MPC问题(20),不妨取Ne=1,δ=1;选取不确定性的边界函数为式(24)描述的单调递减函数,取λ=0.999,ρc=1.4×10-3,后文的分析均选取该边界函数以及当前λ,ρc的值.使用MATLAB实现上述两种算法,其中MPC问题的求解使用非线性最优化函数fmincon,并选用其默认优化算法计算,迭代终止条件取默认值(本文认为当控制和辨识趋于稳定时只需很少的迭代步数就可以获得最优解).所有仿真均在配置为1.7 GHz,Intel Core i5 的处理器的计算机上进行.在上述对比仿真实验中,确定性等价自适应MPC以及本文算法每一步迭代控制的CPU时间分别约为48.27 ms和83.03 ms;实验所用时间数量级相同且,时间差可接受,故本文算法在运行时间上具有可操作性.选取前50步自适应控制的仿真结果绘制图像,如图3所示,其中虚线代表确定性等价自适应MPC算法,实线代表本文具有可参数化不确定性约束的对偶自适应MPC 算法,t表示迭代步数.图3 由式(20)描述的对偶自适应MPC和由式(12)描述的确定性等价自适应MPC,在其他条件相同时的仿真结果对比(其中粗虚线为确定性等价自适应MPC的数据,粗实线为本文算法获得的仿真数据)Fig.3 Simulation results of dual adaptive MPC(solid line)defined by(20)and the CE adaptive MPC(dotted line)defined by(12)图3清楚地体现了本文方法与传统方法的差异.仿真结果表明,两个方法都能随时间推移,使控制输入、输出趋近目标函数中的要求;当辨识的参数值趋于稳定时,两个控制器也很快地在t=8和t=20左右两个时刻分别实现了收敛.对比而言,本文算法在控制的过程中,具有更加平滑的控制输入和控制输出曲线,可以在降低控制损耗的同时提升过程中整体的控制性能.解得的控制输入在一开始距目标零控制输入较远,这是由参数不准及系统不确定性约束引起的算法的正常反应,这种反应使系统具有良好的快速性.本文算法在辨识上的表现与传统方法有一定差异:传统自适应MPC方法辨识较慢收敛,但估计的参数更接近模型真实参数;本文方法参数估计的收敛值距真实值较远,但较快趋于稳定,这加快了控制进程的收敛.参数估计的收敛值远离真实值,实际上是因为更充分地考虑了系统中的不可建模动态;无论辨识的参数有多准确,未建模动态总会对对象真实输出造成影响,使由模型得到的预测输出偏离真实值,所以从另一个角度来说,面向控制的自适应MPC,无需太过准确的参数估计,而应考虑系统中的不确定因素.根据上述分析,可以认为,仿真结果体现了本文具有可参数化不确定性系统的对偶自适应MPC算法的思想和特点:在自适应控制的过程中考虑了系统中存在的不可建模动态,减小了未来的不确定性,从而提升控制的性能.仿真数据证明:本文算法成功。
2021年第5期(总第49卷第363期) No. 5 in 2021 (Total Vol. 49,No. 363 >建筑节能(中英文)Journal of BEE■暖通空调HV&ACdoi : 10.3969/j. issn.2096-9422.2021.05.010基于APC智能控制技术在医院中央空调节能中的应用付磊、罗淋俊2,刘浩2,蔡跃峰1,粘培坤、褚丹雷2(1.厦门大学附属心血管病医院,厦门361006;2.厦门奥普拓自控科技有限公司,厦门361026)摘要:提出了基于高端控制优化算法的中央空调智能化整体解决方案,并详细介绍了方案的系统建 模、控制器设计过程,通过厦门大学附属心血管病医院为研究对象,介绍了该方案的实现与应用,通过验证证明该方案抗干扰性强、控制精度高、非线性处理能力强、鲁棒性能好,能够在提升中央空调系统整体智能化水平的同时实现节能减排。
该技术为中央空调智能化运营提供了一个新的途径。
随着工业互联网技术与数据算法的不断更新,基于A PC多变量控制优化技术的中央空调整体解决方案将得到更为广泛的工程应用。
关键词:先进过程控制;多变量控制;模型预测控制;智能解决方案;实时优化中图分类号:TU83 文献标志码:A文章编号:2096-9422(2021) 05 ■005447Intelligent Solution Scheme of Central Air-conditioning System Based on APCControl Technique in HospitalFU Lei' , LUO Lin-jun2, LIU Hao2, CM Yue-feng', NIAN Pei-kun, CHU Dan-lei2(1. Xiamen Cardiovascular Hospital of Xiamen University,Xiamen 361006, Fujian,China;2.Optimal Process Control Technologies Co.,Ltd.,Xiamen 361026, Fujian,China)Abstract :This paper proposes an intelligent total solution scheme to central air-conditioning based on advanced process control optimization algorithms, introduces the system modelling of the scheme and the design process of controller in detail, and demonstrates the realization and application of the scheme via researching Xiamen Cardiovascular Hospital of Xiamen University. Through verification, the scheme achieves excellent anti-interference performance, high control accuracy, high non-linear processing capacity, and outstanding robust performance, capable of enhancing the overall intelligentization of central air-conditioning system while accomplishing energy conservation and emission reduction. This technology provides a new way for intelligent operation of central air conditioning. Intelligent solution based on APC control technique will be more widely used in engineering applications.Keywords :advanced process control;multi-variable control technique;model predictive control;intelligence solution scheme;real-time optimization〇引言中央空调系统作为近代建筑物温湿度调节的核 心运维系统之一,被广泛应用于星级酒店、大型购物 中心、医院学校、博物馆、无尘车间等,随着人们对舒 适型要求的提高,特别是医院中央空调运行时间长,能耗较大,高峰时占医院总能耗的soy^eo%11~6]。
DMC算法在CSTR温度控制中的应用何美霞;周箩鱼【摘要】连续搅拌反应釜(CSTR)的温度具有强非线性的动态特性,传统的控制方法效果往往不尽人意,而动态矩阵控制 (DMC)算法在处理非线性问题时有一定的优势.以非线性 CSTR 系统的反应釜温度为控制目标,通过仿真试验研究了DMC算法在该系统分别处于理想状态、存在输出干扰及模型失配3种情况时的目标控制效果.仿真结果表明,DMC算法能有效地控制反应釜温度;理想状态下的控制效果表现出上升速度快、调节时间短、无稳态误差等优点;存在输出干扰和系统模型失配时仍有良好的控制效果,具有较强的抗干扰能力和鲁棒性.【期刊名称】《长江大学学报(自然版)理工卷》【年(卷),期】2018(015)005【总页数】5页(P41-45)【关键词】连续搅拌反应釜;温度控制;动态矩阵控制;干扰;模型失配;鲁棒性【作者】何美霞;周箩鱼【作者单位】长江大学电子信息学院,湖北荆州434023;长江大学电子信息学院,湖北荆州434023【正文语种】中文【中图分类】TP273连续搅拌反应釜(Continuous Stirring Tank Reactor,CSTR)是工业过程中广泛使用的一类反应器[1]。
反应器温度及反应物浓度等对产品质量及生产安全有重大影响,因此对这些指标的控制是工业过程控制领域的研究热点。
然而,由于CSTR具有强非线性的动态特性,传统的控制方法(如PID控制、比值控制)效果不尽人意,因此许多学者开始寻求更优的控制方法。
李述清等[2]针对CSTR系统控制问题,设计了一种基于闭环增益成型算法的PID控制器,提高了PID控制器设计的简洁性和鲁棒性。
刘士荣和俞金寿[3]采用神经模糊逆模控制与PID反馈控制相结合的复合控制策略,应用于CSTR的反应物浓度控制中。
通过仿真证明了这类控制策略的有效性和实用性。
刘松等[4]针对CSTR模型设计了具有高增益观测器的非线性鲁棒控制器(ONRC),提出一种简单的控制器的参数整定方法。
模型预测控制设计报告引言模型预测控制(Model Predictive Control,简称MPC)是一种先进的控制算法,它在过程中基于数学模型进行预测,并优化控制动作以使系统的响应最佳化。
本报告将对MPC算法进行介绍,并探讨其在工业控制领域的应用。
MPC算法原理MPC算法的核心思想是通过建立系统的动态模型,预测系统未来的响应,并通过求解优化问题来计算最佳控制动作。
MPC通常包含以下几个步骤:1. 建立数学模型:根据系统的物理特性、实验数据等,建立系统的动态模型。
动态模型可以是线性或非线性的,用差分方程、微分方程、状态方程等形式表示。
2. 预测系统响应:基于已知的系统初始状态和当前的控制动作,利用数学模型进行系统的状态预测。
预测的时间范围可以根据需求进行选择。
3. 优化问题求解:根据预测的系统响应和控制要求,构建一个优化问题,并通过求解优化算法找到最佳的控制动作。
优化问题的目标可以是最小化误差、最大化系统指标等。
4. 调整控制动作:根据求解得到的最佳控制动作,对系统进行调整。
通常需要考虑控制动作的可行性和实时性。
MPC在工业控制中的应用MPC算法在许多工业控制领域中都得到了广泛的应用,并取得了良好的效果。
以下是几个主要的应用领域:1. 化工过程控制:MPC在化工过程控制中的应用十分广泛。
通过准确的模型化和优化求解,MPC能够更好地控制化工过程的温度、压力、浓度等参数,提高产品质量和生产效率。
2. 电力系统控制:MPC在电力系统的控制中也起到了重要的作用。
通过对发电机组的控制,MPC能够减少能量损失、优化电网稳定性,并满足不同的负荷需求。
3. 汽车控制:MPC在汽车控制中被广泛应用于敏感系统(如刹车、悬挂)的控制中。
通过对车辆动力系统的控制,MPC能够提高车辆的操纵性和安全性。
4. 机器人控制:MPC在机器人控制中的应用也逐渐增多。
通过准确的模型预测和动作优化,MPC能够实现机器人的精确控制和路径规划。
第39卷第3期自动化学报Vol.39,No.3 2013年3月ACTA AUTOMATICA SINICA March,2013模型预测控制—现状与挑战席裕庚1,2李德伟1,2林姝1,2摘要30多年来,模型预测控制(Model predictive control,MPC)的理论和技术得到了长足的发展,但面对经济社会迅速发展对约束优化控制提出的不断增长的要求,现有的模型预测控制理论和技术仍面临着巨大挑战.本文简要回顾了预测控制理论和工业应用的发展,分析了现有理论和技术所存在的局限性,提出需要加强预测控制的科学性、有效性、易用性和非线性研究.文中简要综述了近年来预测控制研究和应用领域发展的新动向,并指出了研究大系统、快速系统、低成本系统和非线性系统的预测控制对进一步发展预测控制理论和拓宽其应用范围的意义.关键词模型预测控制,约束控制,大系统,非线性系统引用格式席裕庚,李德伟,林姝.模型预测控制—现状与挑战.自动化学报,2013,39(3):222−236DOI10.3724/SP.J.1004.2013.00222Model Predictive Control—Status and ChallengesXI Yu-Geng1,2LI De-Wei1,2LIN Shu1,2Abstract Since last30years the theory and technology of model predictive control(MPC)have been developed rapidly. However,facing to the increasing requirements on the constrained optimization control arising from the rapid development of economy and society,the current MPC theory and technology are still faced with great challenges.In this paper,the development of MPC theory and industrial applications is briefly reviewed and the limitations of current MPC theory and technology are analyzed.The necessity to strengthen the MPC research around scientificity,effectiveness,applicability and nonlinearity is pointed out.We briefly summarize recent developments and new trends in the area of MPC theoretical study and applications,and point out that to study the MPC for large scale systems,fast systems,low cost systems and nonlinear systems,will be significant for further development of MPC theory and broadening MPC applicationfields. Key words Model predictive control(MPC),constrained control,large scale system,nonlinear systemsCitation Yu-Geng Xi,De-Wei Li,Shu Lin.Model predictive control—status and challenges.Acta Automatica Sinica, 2013,39(3):222−236模型预测控制(Model predictive control, MPC)从上世纪70年代问世以来,已经从最初在工业过程中应用的启发式控制算法发展成为一个具有丰富理论和实践内容的新的学科分支[1−3].预测控制针对的是有优化需求的控制问题,30多年来预测控制在复杂工业过程中所取得的成功,已充分显现出其处理复杂约束优化控制问题的巨大潜力.进入本世纪以来,随着科学技术的进步和人类社会的发展,人们对控制提出了越来越高的要求,不收稿日期2012-06-25录用日期2012-09-29Manuscript received June25,2012;accepted September29, 2012本文为黄琳院士约稿Recommended by Academician HUANG Lin国家自然科学基金(60934007,61074060,61104078)资助Supported by National Natural Science Foundation of China (60934007,61074060,61104078)1.上海交通大学自动化系上海2002402.系统控制与信息处理教育部重点实验室(上海交通大学)上海2002401.Department of Automation,Shanghai Jiao Tong University, Shanghai2002402.Key Laboratory of System Control and Information Processing of Ministry of Education(Shanghai Jiao Tong University),Shanghai200240该文的英文版同时发表在Acta Automatica Sinica,vol.39,no.3, pp.222−236,2013.再满足于传统的镇定设计,而希望控制系统能通过优化获得更好的性能.但在同时,优化受到了更多因素的制约,除了传统执行机构等物理条件的约束外,还要考虑各种工艺性、安全性、经济性(质量、能耗等)和社会性(环保、城市治理等)指标的约束,这两方面的因素对复杂系统的约束优化控制提出了新的挑战.近年来,在先进制造、能源、环境、航天航空、医疗等许多领域中,都出现了不少用预测控制解决约束优化控制问题的报道,如半导体生产的供应链管理[4]、材料制造中的高压复合加工[5]、建筑物节能控制[6]、城市污水处理[7]、飞行控制[8]、卫星姿态控制[9]、糖尿病人血糖控制[10]等,这与上世纪预测控制主要应用于工业过程领域形成了鲜明对照,反映了人们对预测控制这种先进控制技术的期望.本文将在分析现有成熟的模型预测控制理论和工业预测控制技术的基础上,指出存在的问题,综述当前针对这些问题的研究动向,并对模型预测控制未来可能的研究提出若干看法.3期席裕庚等:模型预测控制—现状与挑战2231现有预测控制理论和应用技术存在的问题上世纪70年代从工业过程领域发展起来的预测控制,是在优化控制框架下处理约束系统控制问题的,反映了约束控制的研究从反馈镇定向系统优化的发展.大量的预测控制权威性文献都无一例外地指出,预测控制最大的吸引力在于它具有显式处理约束的能力[1−3,11−12],这种能力来自其基于模型对系统未来动态行为的预测,通过把约束加到未来的输入、输出或状态变量上,可以把约束显式表示在一个在线求解的二次规划或非线性规划问题中.随着预测控制工业应用的普及和软件产品的成熟,标准二次规划算法和序贯二次规划算法被引入预测控制的优化求解.在全球数千个大型工业设施上的成功应用,表明预测控制作为一种实际可用的约束控制算法,已受到了工业过程控制领域的广泛认同[1].Qin等在2003年发表的著名论文[1]中对工业预测控制的发展历程和应用现状做了完整的综述,根据到1999年对于国际上5家主要预测控制软件厂商产品应用的不完全统计,预测控制技术已在全球4600多个装置和过程中得到应用,涉及炼油、石化、化工、聚合、制汽、制浆与造纸等工业领域,预测控制软件产品也已经历了四个发展升级阶段.在我国,预测控制软件开发及典型工程应用被纳入国家“九五”科技攻关,浙江大学、清华大学、上海交通大学等单位都开发了具有自主知识产权的多变量预测控制软件并在一些工业过程中得到成功应用.浙大中控技术有限公司等还实现了预测控制软件的商品化并在国内推广,有力地推动了预测控制在我国的工业应用.尽管预测控制在国内外工业过程中都得到了成功应用,但作为要解决当前经济社会面临的约束优化控制问题的有效技术,仍有以下局限性:1)从现有算法来看,主要还只适用于慢动态过程和具有高性能计算机的环境,从而大大限制了其在更广阔应用领域和应用场合的推广现有的工业预测控制算法需要在线求解把模型和约束嵌入在内的优化问题,每一步都需采用标准规划算法进行迭代,涉及很大的计算量和计算时间,使其只能用于可取较大采样周期的动态变化慢的过程,并且不能应用在计算设备配置较低的应用场合(如DCS的底层控制).Qin在文献[1]中对已投运的线性预测控制产品的应用领域进行了分类,在所统计的2942个案例中,炼油、石化、化工领域占了绝大部分,分别为1985、550、144例.虽然这只是到1999年为止的数据,而且统计的只是国际上主要预测控制商用软件产品的应用状况,但还是趋势性地反映出预测控制的规模应用主要局限在过程工业领域,特别是炼油、石化工业.对于制造、机电、航空等领域内的大量快速动态系统,如果不采用性能较高的计算设备,这类标准优化算法就很难满足小采样周期下的实时计算要求,所以至今未能在这些领域内形成规模应用.2)从应用对象来看,主要还限于线性或准线性过程现有工业预测控制技术的主流是针对线性系统的,成熟的商用软件及成功案例的报道以线性系统为多,虽然软件厂商也推出了一些非线性预测控制产品,但据文献[1]统计,其投运案例数大致只及线性预测控制产品的2%,远未形成规模.即使在过程工业中,预测控制技术的应用也只局限在某些过程非线性不严重的行业,如精炼、石化等,而在非线性较强的聚合、制气、制浆与造纸等领域应用不多.造成这一现象主要是由于在工业过程中非线性机理建模要耗费很大代价,而且很难得到准确的模型,此外非线性约束优化问题的在线求解尚缺乏实时性高的有效数值算法.面对着经济社会发展各行各业对预测控制技术的需求,对象或问题的非线性将更为突出.控制界和工业界都认识到发展非线性预测控制的重要性,例如以非线性模型预测控制为主题的两次国际研讨会NMPC05、NMPC08,就汇聚了国际知名学者和工业界专家认真评价和讨论非线性模型预测控制的现状、未来方向和未解决的问题[13−14].但到目前为止,虽然非线性模型预测控制已成为学术界研究的热点,但在工业实践中仍然处于刚起步的状态[15].3)从应用技巧来看,主要还需依靠经验和基于专用技巧(Ad-hoc)的设计现有的预测控制算法多数采用工业界易于获得的阶跃响应或脉冲响应这类非参数模型,并通过在线求解约束优化问题实现优化控制,对于约束系统无法得到解的解析表达式,这给用传统定量分析方法探求设计参数与系统性能的关系带来了本质的困难,使得这些算法中的大量设计参数仍需人为设定并通过大量仿真进行后验,因此除了需要花费较大的前期成本外,现场技术人员的经验对应用的成败也起着关键的作用,实施和维护预测控制技术所需要的高水平专门知识成为进一步应用预测控制的障碍.30多年来,工业预测控制的技术和产品仍保持着其原有的模式,并没有从预测控制丰富的理论成果中获取有效的支持.最近,应用界已认识到长期以来在过程工业中成功应用但其基本模式保持不变的工业预测控制算法的局限性,研发预测控制技术的著名软件公司Aspen Technology正在考虑摆脱传224自动化学报39卷统的模式,通过吸取理论研究的成果研发预测控制的新产品[16].综上所述,预测控制技术的应用虽然取得了很大的成功,特别在过程控制界已被认为是唯一能以系统和直观的方式处理多变量约束系统在线优化控制的先进技术,但它的应用领域和对象仍因现有算法存在的瓶颈而受到局限,对于更广泛的应用领域和更复杂的应用对象,只能从原理推广的意义上去研究开发相应的预测控制技术,远没有形成系统的方法和技术.此外,现有的工业预测控制算法与近年来迅速发展的预测控制理论几乎没有联系,也没有从中汲取相关的成果来指导算法的改进.因此在解决由于科学技术和经济社会发展所带来的各类新问题时,还面临着一系列新的挑战.与预测控制的实际应用相比较,预测控制的理论研究从一开始就落后于其实践.纵观预测控制理论研究的进程,不难发现它经历了两个阶段[17]:上世纪80年代到90年代以分析工业预测控制算法性能为特征的预测控制定量分析理论,以及上世纪90年代以来从保证系统性能出发设计预测控制器的预测控制定性综合理论.由于后者能够处理包括线性或非线性的对象,包括输入、输出和状态约束在内的相当一般的约束,包括稳定性、优化性能和鲁棒稳定等不同要求的问题,因此引起了学术界极大的兴趣.十多年来,在国际主流学术刊物上已涌现了大量相关论文,呈现出学术的深刻性和方法的创新性,也为约束系统优化控制的研究带来了新的亮点.经过十多年的发展后,预测控制的定性综合理论虽然已取得了丰硕的成果,发表了数以百计的具有很高理论价值的论文,但就目前的研究成果来看,还未能被应用领域所接受.除了这些理论所综合出的算法具有工业界不常采用的模型外,其从综合出发的研究思路也存在着本质的不足.1)物理意义不明确,难以与应用实践相联系预测控制的定性综合理论与定量分析理论不同,在每一时刻的滚动优化中,不是面对一个已有的、根据实际优化要求和约束条件确定的在线优化问题,而需要把在线优化的内容结合控制律一并综合设计.为了得到系统性能的理论保证,往往需要在具有物理意义的原始优化问题中修改性能指标(加入终端惩罚项),加入诸如终端状态约束、终端集约束等人为约束[18],这不但增加了设计的保守性,而且因为这些人为约束与系统受到的实际物理约束一并表达为同一优化问题中的约束条件,使得优化问题中具有物理意义的原始约束湮没在一系列复杂数学公式所表达的整体条件中,这些条件需要通过计算后验,缺乏对实际应用中关注的带有物理意义的分析结论.最典型的如在实际应用中的可行解指的是系统满足所有硬约束的解,而在预测控制定性综合理论中,可行性是指除了满足对系统状态和输入的硬约束外,还要满足包括不变集、Lyapunov函数递减、性能指标上限等在内的由系统设计所引起的一系列附加约束,甚至后者还成为约束的主体,因此很难与应用实践紧密联系.此外,约束下系统状态的可行域有多大,线性矩阵不等式是否有解,如果无解,约束放松到何种程度可以求解等,都无法从现有的研究结果中得到.2)在线计算量大,无法为应用领域所接受预测控制定性综合理论研究的出发点是如何在理论上保证闭环系统在算法滚动实施时的稳定性、最优性和鲁棒性,通常要把原优化问题转化为由新的性能指标和一系列线性矩阵不等式(Linear matrix inequality,LMI)约束描述的优化算法,所以几乎每一篇论文都会根据所研究的问题提出一个甚至多个预测控制综合算法.但是这些研究的重点几乎都放在算法条件如何保证性能的理论证明上,至于算法的具体实施,则认为已有相应的求解软件包即可,并不关注其在线实现的代价.大量人为约束的加入,虽然对系统性能保证是必要的,但同时也极大地增加了优化求解的计算量.特别对鲁棒预测控制问题,由于所附加的LMI条件不但与优化时域相关,而且与系统不确定性随时域延伸的各种可能性有关,LMI的数目将会急剧增长,对在线计算量的影响更为突出.虽然近年来这一问题已开始得到重视,但总体上因其在线计算量大的不足,很难受到应用领域的关注,也很少有在实际中成功应用的案例报道.在预测控制形成的初期,人们曾多次指出其理论研究落后于实际应用,两者之间存在着较大的差距.经过十多年来学术界的努力,虽然形成了成果丰富的预测控制定性综合理论,但由于两者的出发点不同,其理论意义明显高于实用价值,实际上并没有缩小预测控制理论和应用间的差距,远未成为可支持实际应用的约束优化控制的系统理论.综合以上对预测控制应用状况和理论发展的分析可以看出,虽然预测控制的工业应用十分成功,预测控制的理论研究体系也相当完善,但现有的预测控制理论和应用之间存在着严重的脱节,不能满足当前经济社会发展对约束优化控制的要求.我们可以把现有预测控制理论和应用技术存在的问题主要归结为:1)有效性问题.无论是工业预测控制算法还是由预测控制定性综合理论所设计的控制算法,均面临着在线求解约束优化问题计算量大这一瓶颈,极3期席裕庚等:模型预测控制—现状与挑战225大地限制了其应用范围和应用场合.2)科学性问题.预测控制理论研究和实际应用仍有较大距离,商品化应用软件很少吸收理论研究的新成果,理论研究的进展也不注意为实际应用提供指导,缺少既有性能保证又兼顾计算量和物理直观性的综合设计理论和算法.3)易用性问题.目前的预测控制算法都建立在约束优化控制问题一般描述和求解的基础上,对计算环境的要求和培训维护成本都比较高,缺少像PID控制器那样形式简洁、可应用于低配置计算环境、易于理解和掌握的“低成本”约束预测控制器.4)非线性问题.目前预测控制理论和算法的主要成果是针对线性系统的,由于实际应用领域中存在大量非线性控制问题,这方面的研究特别是应用还很不成熟.2当前研究动态随着本世纪科技、经济和社会的发展,各应用领域对约束优化控制的需求日益增长,人们对上面提到的工业预测控制算法和现有预测控制理论的不足有了越来越清晰的认识,促使预测控制理论和应用的研究向着更深的层次发展.当前,模型预测控制已成为控制界高度关注的热点,在各类学术刊物和会议上发表的与预测控制相关的论文居高不下.仅在2007年∼2011年的五年中,通过对Elsevier出版物及IEEE数据库的不完全检索,已查到预测控制相关论文1319篇,其中在Au-tomatica、Control Engineering Practice、Journal of Process Control、IEEE Transactions on Auto-matic Control等刊物上发表的相关论文数分别为74篇、75篇、164篇、35篇.2008年和2011年两次IFAC世界大会上,与预测控制有关的论文分别为131篇和138篇.对预测控制工业应用技术做出全面综述的论文“A survey of industrial model predictive control technology”[1]在2008年IFAC 世界大会上获得CEP最佳论文奖,全面综述预测控制稳定性理论的论文“Constrained model predic-tive control:stability and optimality”[2]在2011年IFAC世界大会上获得了最有影响力奖(High Impact Award).在国内,除了与国际同步开展的对预测控制理论的研究外[19−26],预测控制的应用已从传统的炼油、石化、化工行业延伸到电力[27]、钢铁[28]、船舶[29]、空天[30]、机电[31]、城市交通[32]、渠道[33]、农业温室[34]等领域,各种新的改进算法和策略也屡见报道.通过对近年来国内外预测控制研究工作的分析,可以清楚地看到,一方面,人们对预测控制解决在线约束优化控制寄予很高的期望,试图利用它解决各自领域中更多更复杂的问题;另一方面,工业预测控制算法的不足和现有预测控制理论的局限,又使人们在解决这些问题时不能简单地应用已有的理论或算法,必须研究克服其不足的新思路和新方法.这种需求和现状的矛盾,构成了近年来预测控制理论和算法发展的强大动力,同时也是预测控制理论和算法尽管似乎已很成熟,但人们仍然还在不断研究的主要原因.针对上述预测控制理论和算法的不足,近阶段国内外开展的研究可大致归结为以下几个方面:1)研究降低预测控制在线优化计算量的结构、策略和算法预测控制在线求解约束优化问题计算量大这一瓶颈,极大地限制了其应用范围和应用场合.针对这一问题,人们从结构、策略、算法层面开展了广泛的研究.a)结构层面:递阶和分布式控制结构随着制造、能源、环境、交通、城市建设等领域的迅猛发展,企业集成优化系统、交通控制系统、排水系统、污水处理系统、灌溉系统等大规模系统的预测控制受到了格外的关注[7,35−38],这类大系统的特点是组成单元多、模型复杂、变量数目巨大,整体求解其大规模约束优化问题在实际中几乎不可行.因此,针对实际系统的应用需求,人们普遍借鉴传统大系统理论提供的递阶控制结构把整体优化求解的复杂性进行分解.虽然基于同一模型分解协调的多级递阶控制方法在理论上已发展得较为成熟,但考虑到模型和实际环境的复杂性,在研究中通常更倾向于应用在不同层次采用不同模型的多层递阶结构[39],其研究的重点在于确定各层次的模型和优化控制目标以及协调各层次之间的关系,由此发展有效可行的控制框架和算法,所提出的控制方案和算法常通过仿真或实际运行数据加以验证.在大规模系统预测控制的研究中,近年来更受到重视的是采用分布式结构降低计算复杂性[40−41],分布式预测控制基于用局部信息进行局部控制的思想把大规模约束优化控制问题分解为多个小规模问题,不仅可以大大降低计算负担,而且提高了整体系统的鲁棒性.分布式预测控制的研究重点包括各子系统之间耦合关联的处理、子系统的优化决策及相互间的信息交换机制、全局稳定性的保证及最优性的评估等[42].近年来通信技术的发展和分布式控制软硬件的完善,使分布式预测控制从理论走向实践,应用已遍及到多个领域,包括过程控制[43]、电力系统[44]、交通系统[45]及近年来十分活跃的多智能体协作控制等[46].b)策略层面:离线设计/在线综合与输入参数226自动化学报39卷化策略在预测控制定性综合理论研究中,虽然系统性能可得到严格的理论保证,但设计所带来的额外计算负担十分庞大,导致本来已成为应用瓶颈的在线计算复杂性更为突出,这也是应用界对这些理论研究成果可用性的主要质疑.针对这一问题,在预测控制的定性综合中提出了“离线设计、在线综合”的策略,通过把所综合控制律的部分在线计算转换为离线计算,达到降低在线计算量的目的.文献[47]应用该思路给出了文献[18]提出的约束鲁棒预测控制器的简化设计方法;文献[48]利用名义系统指标离线设计不变集序列,在线时通过核算当前状态所在的最优不变集来确定控制律;类似的设计还包括文献[49];文献[50]通过离线求解有限时域优化控制序列,并采用Set membership来得到近似最优解,以提高求解效率.在这里特别要提到的是由Bemporad等提出的显式(Explicit)模型预测控制器[51−52],它通过对预测控制在线约束优化问题的分析,离线求解多参数规划问题,对约束状态空间分区并设计各区间的显式反馈控制律;在线控制时,只需依据系统的当前状态,选择实施相应的状态反馈控制律.这种方法把大量计算转移到离线进行,在线控制律的计算十分简易,而且有坚实的理论基础,因此受到了广泛的关注,进一步研究算法简化和对非线性系统的推广、以及算法在微处理器中的应用等也已见报道,如文献[53−54].但该方法离线需求解一个NP-hard的多参数规划问题,离线计算量随着问题规模增大而急剧增加,同时由于分区数的指数增长而导致巨大的内存需求,只能应用于小规模的问题[55].为此,近年来国内外学者进行了进一步的探索.文献[56]采用分段连续网格函数(Lattice PWA function)表示显式预测控制的解,以降低其对于存储空间和在线计算的要求;文献[57]通过分析二次规划问题求解方法在存储和计算方面的复杂度,提出一种半在线的显式预测控制算法,在存储量和在线计算时间之间进行平衡;文献[58]将动态规划和显式预测控制方法相结合,把预测控制的优化问题分解为小规模问题;而文献[59]针对非线性系统预测控制问题提出了近似的显式预测控制方法.离线设计、在线综合的方法能有效地解决预测控制在线优化计算量大的瓶颈问题,但要求原有的预测控制器设计方法可以进行分解,并且需要为在线综合保留一定的自由度,因此不能适用于所有的预测控制定性综合算法.在工业预测控制算法中,为了降低在线优化的计算量,很早就采用了启发式的“输入参数化”策略[1],包括输入“分块化(Blocking)”技术[60]和预测函数控制算法[61],前者把一定时间段内的控制量设置为不变,以减少控制自由度的代价来降低在线优化问题的规模,后者则把控制量表达为一组基函数的组合,使在线优化变量转化为数目较少的基函数的系数.这些策略虽然有很强的实用性并已大量应用于实际过程,但缺乏对系统稳定性等的理论保证.现有的预测控制稳定性综合方法在用于这类系统时,又因输入参数化造成递归可行性难以保证而不能奏效.近年来国内外学者对此进行了进一步的研究.对于Blocking技术,文献[62]采用时变的集结矩阵保证集结预测控制器的闭环稳定性,文献[63−64]就集结预测控制器的可行性问题进行了研究,并提出改善其可行性的方法.文献[65]提出了预测控制优化变量的广义集结策略,这一框架不但可以涵盖以上两种方法,而且由于把原有输入参数化的物理映射扩展为集结变换的数学映射,提供了更大的设计自由度,也为系统分析建立了必要的基础.在此基础上,文献[66]进一步设计了等效/拟等效集结策略以改善集结预测控制的控制性能.c)算法层面:各种改进或近似优化算法针对约束预测控制在线优化的问题形式,对标准优化算法进行改进或做适当近似,也是近期来降低预测控制在线计算量的一类尝试.文献[3]列举了在线求解大型二次规划(Quadratic programming, QP)和非线性规划问题时对算法的若干改进工作.文献[55]提出用扩展的Newton-Raphson算法取代现有算法中常用的二次规划和半定规划(Semi-definite programming,SDP)算法,可使计算量降低10倍以上.文献[67]提出了三种针对预测控制在线求解QP问题的快速算法.文献[68]打破了求解优化问题中“优化直至收敛”的概念,提出了实时迭代的概念,它在每一采样时刻只需计算一次迭代,其结果通过特定的移位与下一时刻的优化问题联系起来,在此基础上,文献[69]又提出了基于伴随导数和非精确雅可比阵的优化算法.文献[70]提出了采用部分列表的快速、大规模模型预测控制方法.此外,采用神经网络求解二次规划等问题又有了新的发展,与以往的工作相比,新的神经网络方法在保证收敛到全局最优解及降低神经网络结构复杂度方面都取得了较好的结果[71].2)鲁棒预测控制理论的研究更加注重实际可用性鲁棒预测控制理论在上世纪90年代中期已初步形成,从本世纪以来更成为预测控制理论研究的重点,在已有大量成果基础上,近年的研究更注意向实际靠拢和解决相关的难点问题.。
Chance-constrained model predictive control for spacecraft rendezvous with disturbance estimationFrancisco Gavilan a,Rafael Vazquez a,n,Eduardo F.Camacho ba Departamento de Ingenierı´a Aeroespacial,Escuela Superior de Ingenieros,Universidad de Sevilla,Camino de los Descubrimientos s/n,41092Sevilla,Spainb Departamento de Ingenierı´a de Sistemas y Automa´tica,Escuela Superior de Ingenieros,Universidad de Sevilla,Camino de los Descubrimientos s/n,41092Sevilla,Spaina r t i c l e i n f oArticle history:Received17December2010Accepted23September2011Keywords:Predictive controlRobust controlBounded disturbancesConstraint satisfaction problemsSpacecraft autonomySpace vehiclesa b s t r a c tA robust Model Predictive Controller(MPC)is used to solve the problem of spacecraft rendezvous,usingthe Hill–Clohessy–Wiltshire model with additive disturbances and line-of-sight constraints.Since astandard(non-robust)MPC is not able to cope with disturbances,a robust MPC is designed using achance-constrained approach for robust satisfaction of constraints in a probabilistic sense.Disturbancesare modeled as Gaussian allowing for an explicit transformation of the probabilistic constraints intosimple algebraic constraints.To estimate the distribution parameters a predictor of disturbances isproposed.Both robust and non-robust MPC control laws are compared using the Monte Carlo method,which shows the superiority of the robust MPC.&2011Elsevier Ltd.All rights reserved.1.IntroductionTechnology enabling simple autonomous spacecraft rendezvousand docking is becoming a growingfield as access to spacecontinues to increase.After decades of development,manyapproaches have been proposed and there have been many experi-ences,positive and negative;see Woffinden and Geller(2008)for anhistorical account or Fehse(2003)for the basics.For instance,one ofthe most recent developments in thefield is ESA’s AutomatedTransfer Vehicle(ATV),mainly developed by EADS Astrium,anexpendable unmanned spacecraft designed to resupply the Interna-tional Space Station.ATV has automatic rendezvous capabilities,asdemonstrated in itsfirst successfulflight in2008.Thefield has become very active in recent years,with anincreasingly growing literature;for instance,among many,onecan cite Richards,Schouwenaars,How,and Feron(2002),wherefuel-optimal trajectories with avoidance constraints are designedusing mixed-integer linear programming,Wang,Mokuno,andHadaegh(2003),which includes autonomous rendezvous anddocking capabilities into formationflying satellites,Geller(2006),which uses a linear covariance analysis method to designimpulsive maneuvers,or Breger and How(2008),where safe,fail-tolerant rendezvous trajectories are planned.This work approaches the problem of rendezvous of spacecraftusing a chance-constrained Model Predictive Control(MPC)withon-line prediction of disturbance statistical properties.MPC(Camacho&Bordons,2004)originated in the late seven-ties and has developed considerably since then.There are manyapplications of predictive control successfully in use at thecurrent time,not only in the process industry but also in otherapplications ranging from solar technology(Camacho,Berenguel,&Bordons,1994)toflight control(Breger&How,2006).ModelPredictive Control is considered to be a mature technique forlinear and rather slow systems like the ones usually encounteredin the process industry.The term Model Predictive Control does not designate aspecific control strategy but rather an ample range of controlmethods which make explicit use of a model of the process toobtain the control signal by minimizing an objective function overafinite receding horizon.In MPC the process model is used topredict the future plant outputs,based on past and current valuesand on the proposed optimal future control actions.These actionsare calculated by the optimizer taking into account the costfunction(where the fuel cost and the future tracking error areconsidered)as well as the constraints.One of the advantages of MPC is that robust control methods canbe easily incorporated.In space vehicles,one canfind multiplesources of disturbances,such as position or velocity measuring errors,thruster misalignments,or even atmospheric drag;so there is a needto design robust control schemes to deal with these disturbances.Thus,MPC is very suitable to deal with the problem of space-craft rendezvous,which is inherently slow and can be veryprecisely modeled by linear equations(shown in Section2).Theuse of robust MPC for rendezvous of spacecraft is not new;forinstance Richards and How(2003)analyzes the advantages ofrobust and non-robust MPC for rendezvous compared with otherContents lists available at SciVerse ScienceDirectjournal homepage:/locate/conengpracControl Engineering Practice0967-0661/$-see front matter&2011Elsevier Ltd.All rights reserved.doi:10.1016/j.conengprac.2011.09.006n Corresponding author.Fax:þ34954486041.E-mail addresses:fgavilan@us.es(F.Gavilan),rvazquez1@us.es(R.Vazquez),eduardo@.es(E.F.Camacho).Control Engineering Practice](]]]])]]]–]]]methods.In the work of How and Tillerson(2001)the effect of velocity measurements error during formationflight is taken into account.Sensor errors are modeled,and a robust MPC scheme is proposed that satisfies the constraints for the worst case dis-turbance,recalculating the trajectory only when the spacecraft get out of a desired error box.In these works,the key idea is to explicitly take into account system disturbances and uncertainties and to optimize the objective function for the worst case scenario(Camacho& Bordons,2004).However,for these methods it is necessary to obtain an estimate on the bound of the disturbances.In Richards (2004),several methods for estimation of uncertainty properties are proposed.This work proposes the use of the so called chance-constrained model predictive control as an alternative to design robust control-lers for spacecraft.In this approach,disturbances are incorporated into the problem constraints using a probabilistic formulation.A procedure to transform these probabilistic constraints into algebraic equations is given,and control signals are computed so the con-straints are satisfied with a desired probability.Chance-constrained MPC can be found in previous works(mainly for chemical engineer-ing applications),for instance,in Li,Wendt,and Wozny(2000)and Schwarm and Nikolaou(1999)where methods are proposed to deal with linear systems in which uncertainties are present in the step response coefficients of the systems.These authors consider that statistical properties of unknown parameters are acquired off-line.This work employs an on-line estimator of statistical properties of disturbances together with the chance-constrained formulation of the problem.The advantage of the chance-constrained formulation is that it does not need to know a priori bounds on the size of the disturbances.However,it needs a distribution model for them;a Gaussian model is used which allows for an explicit solution of the probabilistic constraints into algebraic constraints,thus allowing for fast computation of a solution.The parameters of the Gaussian model are inferred from past disturbances using the on-line estimator.To the best of our knowledge,these methods have not been proposed before to solve the problem of spacecraft rendezvous.The structure of this work is as follows.Section2describes the mathematical model for rendezvous spacecraft used for MPC and the constraints of the rendezvous problem.Next,Section3follows with a formulation of standard(non-robust)Model Predictive Control suitable for the rendezvous maneuver with continuous thrust.Then the robust chance-constrained MPC is formulated with estimation of disturbance properties.Section4shows a Monte Carlo comparison of the robust and non-robust methods.The comparison is also shown for elliptical target orbits,with the discrepancies due to eccentricity considered as a disturbance.Section5closes the work with somefinal remarks.2.Model of spacecraft rendezvousThere are numerous mathematical models for spacecraft rendezvous;which one should be used depends on the para-meters of the scenario.In Carter(1998)a survey of numerous mathematical models for spacecraft rendezvous can be found.For instance,if the target is orbiting in a circular Keplerian orbit,the general equations of the relative movement between an active chaser spacecraft and a passive target vehicle are (see Wie,1998)€x¼2n_yþn2ðRþxÞÀm Rþx½ðRþxÞþy2þz2 þu x,€y¼À2n_xþn2yÀm y½ðRþxÞ2þy2þz2 3=2þu y,€z¼Àm z½ðRþxÞ2þy2þz2 3=2þu z,ð1Þwhere x,y,and z denote the position of the chaser in a local-vertical/local-horizontal(LVLH)frame of referencefixed on thecenter of gravity of the target vehicle(see Fig.1),in which x refersto the radial position,y to the in-track position,and z to the cross-track position.The velocity of the chaser in the LVLH frame isgiven by_x,_y,and_z;and the variables u x,u y,and u z are the inputs(thrust actuation)acting on the chaser vehicle.R is the target orbitradius and n¼ffiffiffiffiffiffiffiffiffiffiffim=R3qis the angular speed of the target throughits orbit(where m is the gravitation parameter of the Earth,m¼398600:4km3=s2).Moreover,if the approaching vehicle is close to the target,Eq.(1)can be linearized around the rendezvous position,leadingto the linear Hill–Clohessy–Wiltshire(HCW)equations(intro-duced in Hill,1878and Clohessy&Wiltshire,1960)whichdescribe with adequate precision the relative position of thespacecraft.The HCW model is the one used throughout thispaper,including the possibility of disturbances to allow forunmodeled effects,and the rotation of the target vehicle.It must be noted that,in many situations,the HCW equationsare not accurate.For instance,if the target vehicle is moving ina Keplerian eccentric orbit(see Inalhan,Tillerson,&How,2002)or if some orbital perturbations are taken into account(see forexample Humi&Carter,2008).Section4considers simulationswith the target orbiting in an eccentric Keplerian trajectory andshows that the control design(based on the HCW equations)still works.Considering that the control inputs are constant for eachsample time interval of duration T,it is possible to derive thefollowing discrete time version of the HCW equations:x kþ1¼Ax kþBu kþd k,ð2Þwhere an unknown vector d k has been added to take into accountpossible additive disturbances.In(2),x k,u k and d k denote,respectively,the state,input,anddisturbance at time k,wherex¼½x y z_x_y_z T,u¼½u x u y u z T,ð3Þd¼½d x d y d z d_x d_y d_z T,ð4Þwhere d x,d y,d z,d_x,d_y,and d_z represent the disturbances enteringthe system.Both are referred to the LVLH axes as indicated bytheir respectivesubscripts.Fig.1.LVLH frame.F.Gavilan et al./Control Engineering Practice](]]]])]]]–]]]2The system(2)will be used for predicting the spacecraft position in the predictive controller formulation(Section3).The matrices A and B appearing in(2)are given byA¼4À3C00Sn2ð1ÀCÞn6ðSÀnTÞ10À2ð1ÀCÞn4SÀ3nTn00C00Sn3nS00C2S0À6nð1ÀCÞ00À2S4CÀ3000ÀnS00C2666666666437777777775,ð5ÞB¼1ÀCn22nTÀ2Sn22ðSÀnTÞn2À3T22þ41ÀCn2001ÀCn2Sn21ÀCn2ðCÀ1ÞnÀ3Tþ4Sn00Sn26666666666643777777777775,ð6Þwhere S¼sin nT and C¼cos nT.The disturbances are unknown,so it is assumed that d k is a random vector,with known distribution but unknown distribution parameters mean d and covariance R, respectively.These disturbances might arise from errors in the input signals(as thrusters are typically subject to command uncertainties and are never perfectly aligned),or they could also be thought of as unmodeled dynamics(in this case they are not random;however, the randomness assumption is kept for convenience).In Section4.3 the disturbance model used in simulations is described.Even though the disturbances are modeled as additive,in Section4.3it is shown that the control scheme works for other kind of disturbances such as multiplicative disturbance or mod-eling errors.2.1.Constraints on the problemTwo set of constraints are considered in this paper.First,for sensing purposes(see Breger&How,2008)it is required that the chaser vehicle remains inside a line of sight(LOS)area from the docking point;and second,the amount of thrust that can produce the actuators is bounded.The LOS area is a region defined to guarantee that the chaser spacecraft is all time visible from the docking point.Thus this area must be defined using a new bodyfixed frame,since the target can rotate respect to the LVLH axes used in(2),which arefixed to the orbit.Then,once the LOS region is formulated in body axes,a transformation must be used to include these constraints into the rendezvous problem,which is formulated in the LVLH frame.The target bodyfixed reference frame is shown in Fig.2.In this reference system,one can define the LOS region by the equations y B Z c xðx BÀx0Þ,y B ZÀc xðx Bþx0Þ,y B Z c zðz BÀz0Þ,y B ZÀc zðz Bþz0Þand y B Z0(where x B,y B and z B denote the coordinates in the bodyfixed frame);as shown in Fig.2.The LOS constraint is formulated as A L x B k r b L,whereA L¼0À10000c xÀ10000Àc xÀ100000À1c z0000À1Àc z0002666666437777775,ð7Þb L¼½0c x x0c x x0c z z0c z z0 Tð8Þand x B k denotes the state in the bodyfixed reference frame.Since these constraints are not defined using the same refer-ence frame than the equation of motion used by the controller(2), a transformation of LOS constraints from body axes to LVLH frame must be done.The transformation of these matrices can be easily computed using projective geometry(see Hartley&Zisserman, 2003).Using homogeneous coordinates,one can write the equations of a set of n planes aspðnÂ4Þ~x¼0,ð9ÞFig.2.Line of sight region.F.Gavilan et al./Control Engineering Practice](]]]])]]]–]]]3where each row of p ,namely p i ,defines one plane;and ~xis the vector of homogeneous coordinates,that is ~x¼½x y z 1 T :ð10ÞIt can be proven that under a projective transformation ~x¼H ~x 0,a plane transforms as p 0i ¼p i H ,ð11Þwhere H is a transformation matrix in homogeneous coordinates.In this case,the LOS planes in body axes introduced in (7)and (8)can be defined as:0À100c xÀ10Àc x x 0Àc x À10Àc x x 00À1c z Àc z z 00À1Àc z Àc z z 02666666437777775|fflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflffl{zfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflffl}p Bx By B z B 12666437775|fflfflffl{zfflfflffl}~x B¼0:ð12ÞThe projective transformation between body axes and LVLH frame can be defined as~x B ¼R 3Â3t 3Â101Â31 !|fflfflfflfflfflfflfflfflfflfflfflfflffl{zfflfflfflfflfflfflfflfflfflfflfflfflffl}H~x LVLH ,ð13Þwhere R is a rotation matrix from LVLH frame to body axes,and tis a translation vector which contains the coordinates of center of the LVLH frame respect to the center of the body axes.Thus,the set of constraint planes in the LVLH frame can be computed as:p LVLH ¼p B Hð14ÞFor instance,if the target spacecraft is rotating around thez LVLH axis with angular velocity O ,the transformation matrix H is defined asH ¼C O S O 00ÀS O C O 000010012666437775,ð15Þwhere C O ¼cos ðO t Þand S O ¼sin ðO t Þ.Thus,the transformed LOS lines are computed as follows:ÀS O ÀC O 00c x C O ÀS O Àc x S O ÀC O 0Àc x x 0Àc x C O ÀS O c x S O ÀC O 0Àc x x 0ÀS O ÀC O c x Àc x z 0ÀS O ÀC O c z Àc z z 02666666437777775|fflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflffl{zfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflffl}p LVLHx y z 12666437775¼0,ð16Þwhere now x ,y and z are the coordinates of the chaser spacecraftin the LVLH frame.Using these terms,the LOS constraint matrices become:A 0L ¼ÀS O ÀC O 0c x C O ÀS O Àc x S O ÀC O 0Àc x C O ÀS O c x S O ÀC O 0H 5Â3ÀS O ÀC O c z ÀS O ÀC O c z 2666666437777775,b 0L ¼½0c x x 0c x x 0c z z 0c z z 0 T ,ð17Þwhere H 5Â3is a matrix full of zeros,with dimensions 5Â3.Notice that in this situation A 0L and b 0L become time dependent.Then the LOS constraint in the LVLH frame can be rewritten as:A 0L x r b 0L :ð18ÞDealing with the control inputs constraints,it is assumed thatthey are bounded above and below u min r u k r u max ,ð19Þand that u k can take any value in the interval,i.e.,it is assumedthat thruster valves can be opened partially to produce the exact amount of force.3.Robust MPC formulationNext a robust MPC scheme is formulated;first some notation is developed to formulate the general problem,and afterwards it is explained how to tackle the disturbances appearing in (2).3.1.Prediction of the stateThe state at time k þj ,given the state at time k ,and the input signals and disturbances from time k to time k þj À1,is computed by applying recursively Eq.(2):x k þj ¼A jx k þX j À1i ¼0Aj Ài À1Bu k þi þX j À1i ¼0A j Ài À1d k þi :ð20ÞDefine now x S ðk Þ,u S ðk Þ,d S ðk Þas a stack of N p Án x states,N p Án u input signals,and N p Án x disturbances beginning at time k ,where N p is the prediction horizon,n x is the number of state variables and n u is the number of control inputs:x S ðk Þ¼x k þ1x k þ2^x k þN p 266664377775,u S ðk Þ¼u k u k þ1^u k þN p À1266664377775,d S ðk Þ¼d kd k þ1^d k þN p À1266664377775:ð21ÞThen,x S ðk Þ¼Ax k þBu k þd kA 2x k þX 1i ¼0A 1Ài Bu k þi þd k þiÀÁ^A N p x k þX N p À1i ¼0A N p À1Ài ðBu k þi d k þi Þ266666666664377777777775,ð22Þwhich can be written as x S ðk Þ¼Fx k þG u u S ðk ÞþG d d S ,ð23Þwhere G u and G d are block lower triangular matrix with its non-null elements defined by ðG u Þij ¼A i Àj B and ðG d Þij ¼A i Àj (with i Z j ),and the matrix F is defined as:F ¼A A 2^A N p2666437775:ð24ÞNote that it is assumed that the control logic has perfect knowledge of the state vector x k .3.2.Objective functionTaking mathematical expectation,define ^xk þj 9k ¼E ½x k þj ,the expected value of x k þj given x k .Similarly define ^xS ðk þj 9k Þ¼E ½x S ðk þj Þ .For the MPC formulation the following cost functionF.Gavilan et al./Control Engineering Practice ](]]]])]]]–]]]4is used:JðkÞ¼X N pi¼1½^x Tkþi9kRðkþiÞ^x kþi9k þX N pi¼1½u TkþiÀ1Id3Â3u kþiÀ1 ,ð25Þwhere the matrix RðkÞis defined asRðkÞ¼g hðkÀk aÞId3Â3H3Â3 H3Â3H3Â3"#:ð26ÞIn(26),h is the step function,k a is the desired arrival time for docking,g a large positive number,and Id3Â3,H3Â3are,respec-tively,the identity matrix and a matrix full of zeros,both of order3by3.The reason for choosing(26)is that it is desired to arrive at the origin at time k a(and remain there)and at the same time minimize the control effort.Remark1.The relative position during the maneuver is not impor-tant as long as constraints are satisfied and docking is reached on time.Using(23),and since E½dðkþiÞ ¼d,Eq.(25)can be rewritten as:JðkÞ¼ðG u u SðkÞþFx kþG d d SÞT R SðG u u SðkÞþFx kþG d d SÞþu T S Q S u S,ð27Þwhere d S is a stack vector with d repeated N p times,Q S¼Id3NpÂ3N p and R S is a block diagonal matrix defined by:R S¼Rðkþ1Þ&RðkþN pÞ264375:ð28ÞUsing the notation above developed with the LOS constraints formulated in Section2.1,the constraints equations for the state can be rewritten as:A c x S r b c,ð29Þwhere A c and b c are given by:A c¼A0Lðkþ1Þ&A0LðkþN pÞ264375,b c¼½b0L ðkþ1ÞÁÁÁb0LðkþN pÞ T:ð30ÞUsing Eq.(23),one can reformulate the LOS constraints as constraints for the control signals in the following way:A c G u u S r b cÀA c Fx kÀA c G d d S,ð31Þand similarly it is possible to write(19)as:~uminr u S r~u max,ð32Þbeing~u min and~u max stacks of N pÂn u minimum and maximum bounds of the control input.putation of control inputTo compute the control input at time k,one seeks the control signal that minimizes the cost function over the prediction horizon,satisfying at the same time the constraints:minu SJðx k,u S,d SÞ,s:t:A c G u u S r b cÀA c Fx kÀA c G d d S8d S,~uminr u S r~u max:ð33ÞSince the cost function is quadratic and the constraints are linear,if the future disturbance d is perfectly known(for example, in the undisturbed case)then(33)can be solved;the control u k is set to thefirst three components of u S,and the computation is repeated for every time step.However,if the disturbances are not known but rather their mean and variance are known,it is necessary to modify(33),as it is explained next.3.4.Robust satisfaction of constraintsTo eliminate the disturbances from(31),a bound of the term ÀA c G d d S must be found to enforce the satisfaction of constraints in presence of disturbances.Two procedures are given,depending on which disturbance properties are available.Assumefirst that some bounds for the disturbances are known, i.e.,d has the property thatðd xÞmin r d x rðd xÞmax and similarly for the rest of the components of d.Those bounds are summarized as A d d S r c d.Hence,it is assumed that the region defined by this constraint is enclosed by a polytope.Then,it is possible to eliminate the disturbances from(31)by bounding the termÀA c G d d S.This would giveA c G u u S r b cÀA c Fx kÀA c G d d S r b cÀA c Fx kþb dðkÞ,ð34Þwhere b dðkÞis column vector,whose i-th termðb dðkÞÞi is given byðb dðkÞÞi¼mind Sa i d S,s:t:A d d S r c d,ð35Þwhere a i is the i-th row of the matrixÀA c G d.Since the function to minimize is linear and the feasible region is enclosed by a polytope,this minimization can be rapidly solved.Eq.(34)represents the constraints computed for the worst-case disturbances.Hence,enforcing(34)the constraints(29)are robustly satisfied,i.e.,satisfied for any possible disturbance.However,if the disturbance bounds are not precisely known, but the disturbance is modeled as a random vector,the inequality b dðkÞrÀA c G d d S is made to be satisfied with a certain probability. This probability should be near one,thus guaranteeing that the inequality is satisfied for almost all disturbances.Assuming that the disturbances are normally distributed (d$N6ðd,RÞ),and that their mean(d and covariance matrix (R¼R T40)are known,it is possible to write(see Rencher,1998) d$N6ðd,RÞ)ðdÀdÞT RÀ1ðdÀdÞ$w2ð6Þ,ð36Þwhere w2ð6Þis a chi-square probability distribution with six degrees of freedom.Assuming that the statistical properties of the disturbances are time-invariant,Eq.(36)is valid for the disturbances at all times kþj for j¼0,...,N pÀ1:ðd kþjÀd T RÀ1ðd kþjÀd$w2ð6Þ,j¼0,...,N pÀ1:ð37ÞHence,finding a from the equation:Pðw2ð6Þr aÞ¼p,ð38Þit is guaranteed with probability p thatðd kþjÀd T RÀ1ðd kþjÀd r a:ð39ÞThus,p is a parameter of the control design and should be close to unity.Then,dividing the inequality by a,b dðkÞcan be found by solving the following minimization problem for each row i of A c G d.ðb dðkÞÞi¼mind Sa i d S,s:t:ðd kþjÀdÞTða RÞÀ1ðd kþjÀdÞr1,j¼0,...,N pÀ1,ð40ÞF.Gavilan et al./Control Engineering Practice](]]]])]]]–]]]5where a i andðb dðkÞÞi are the i-th row of the matrixÀA c G d and the vector b dðkÞ,respectively.Now,breaking down the stack vector d S into each of itscomponents d k up to d kþNp À1,it is possible to write a i d S¼P N pÀ1j¼0a ij d kþj.Thus,the minimization problem can be written as:ðb dðkÞÞi¼mind kþj XN pÀ1j¼0a ij d kþj,s:t:ðd kþjÀdÞTða RÞÀ1ðd kþjÀdÞr1,j¼0,...,N pÀ1:ð41ÞDefining zðjÞ¼H12ðd kþjÀd,where H¼ða RÞÀ1(being H¼H T40),it is possible to write(41)as:ðb dÞi¼minzðjÞXN pÀ1j¼0ða ij HÀ1=2zðjÞþa ij d,s:t:zðjÞT zðjÞr1,j¼0,...,N pÀ1,ð42Þwhich can be rewritten asðb dÞi¼XN pÀ1j¼0minzðjÞða ij HÀ1=2zðjÞÞþa ij d,s:t:zðjÞT zðjÞr1,j¼0,...,N pÀ1:ð43ÞProblem(43)can be explicitly solved independently for each j via the Lagrange formalism,yielding the minimum atz nðjÞ¼ÀHÀ1=2a Tijffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffia ij HÀ1a Tijq:ð44ÞSubstituting into(43)the rows of the vector b dðkÞareðb dðkÞÞi¼XN pÀ1j¼0Àffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffia ij HÀ1a Tijqþa ij d:ð45ÞOnce the vector b dðkÞis calculated(using Eq.(35)or(45)),the control input at time k is now computed fromminu SJðx k,u S,d SÞ,s:t:A c G u u S r b cÀA c Fx kþb dðkÞ,~uminr u S r~u max,ð46Þwhere now everything is known except for the control inputs to be computed.In simulations,the probabilistic method is used to compute the robust MPC control.3.5.Disturbance estimation algorithmThe robust satisfaction of constraints presented in Section3.4 requires knowledge of disturbance statistical properties.How-ever,it is often the case that such properties are completely unknown and have to be obtained on-line.To do so,it is assumed that the disturbances are normally distributed with mean d and covariance matrix R,i.e.,d$N6ðd RÞ.At each time k,d and R are estimated taking into account all past disturbances,which can be computed a posteriori asd i¼x iþ1ÀAx iÀBu i,ð47Þfor i¼1,...,kÀ1.Then^d k and^R k,the estimates of d and variance R at time k, based on disturbances up to time kÀ1,are given by^d k ¼P kÀ1i¼0eÀlðkÀiÞd iP kÀ1i¼0eÀlðkÀiÞ,ð48Þ^Rk¼P kÀ1i¼0eÀlðkÀiÞðd iÀ^d kÞðd iÀ^d kÞTP kÀ1i¼0eÀlðkÀiÞ,ð49Þwhere l40is a forgetting factor.Even though it has been assumedthat the disturbances are just a random variable,this would helpaccommodate the case in which they are a random process,i.e.,their statistical properties change with time.Define g k¼P kÀ1i¼0eÀlðkÀiÞ.Using the sum of a geometric pro-gression,gk¼eÀlð1ÀeÀl kÞ1ÀeÀl:ð50ÞThen,it is possible to define recursive formulas for(48)and(49)as follows:^dk¼eÀlgkðg kÀ1^d kÀ1þd kÀ1Þ,ð51Þ^Rk¼eÀlgkðg kÀ1^R kÀ1þðd kÀ1À^d kÞðd kÀ1À^d kÞTÞ,ð52Þwith^d0¼0,^R0¼0.These formulas allow to save memory,only needing to storethe last estimate of the mean and covariance.Once the mean and covariance are estimated,the procedureoutlined in Section3.4can be used.4.Simulation results4.1.Rendezvous modelIt is important to remark that although the controller shown inSection3is designed using the linear HCW model(2),in simula-tions the general nonlinear model of spacecraft rendezvous(1)has been considered.Dealing with the model parameters,it has been considered thatthe case of a target spacecraftflying in a circular orbit at500kmof altitude,which means that R0¼6878km and n¼1:1068Â10À3rad=s.Concerning the constraints,the maximum and minimumamount of acceleration that can provide the chaser’s actuatorsare u max¼10À3m=s2and u min¼À10À3m=s2,respectively(all theactuators have the same values).The LOS area is estated with theparameters:x0¼z0¼1:5m and c x¼c z¼1(see Fig.2).−5051015205101520x [m]y[m]XY plane trajectoryFig.3.Non-robust MPC without disturbances(solid line)and with disturbances(dashed line).For clarity,only the XY plane is shown.F.Gavilan et al./Control Engineering Practice](]]]])]]]–]]]6。