1、第七章第七章 最优控制的计算方法最优控制的计算方法本章主要内容7.1 直接法 7.2 间接法 7.3 小结 返回主目录 在前面讨论变分法、极小值原理和动态规划时,我们列举了一些例子。为了易于说明问题,这些例子都是非常简单的,可以用手算来解决问题。但是在实际工作中所遇到的最优控制问题,一般都是很复杂的,必须用计算机求解。因此,最优控制的计算方法就变得十分重要了。这方面的内容十分丰富,由于篇幅所限,我们只介绍几种典型的算法。(无约束)(ii)哈密顿函数 取极小的必要条件H0UHU(有约束)(7-2)或U),(),(min*tUXHtUXHU 由极小值原理可知,最优控制问题的解必须满足以下几个条件(
2、iii)边界条件(包括横截条件)HXXH(i)正则方程(7-1)最优控制的计算方法一般是先求出满足上面三个条件中某两个的解,然后用合适的迭代计算形式逐次改变这个解,以达到满足剩下的另一个条件的解(即最优解)。通常把最优控制的计算方法分成两类:直接法和间接法。p直接法。它的特点是,在每一步迭代中,不一定要满足 取极小的必要条件,而是逐步改善它,在迭代终了使它满足这个必要条件,而且,积分状态方程是从 到 ,积分协态方程是从 到 ,这样就避免了去寻找缺少的协态初值 的困难。常用的直接法有梯度法,二阶梯度法,共轭梯度法。)(tUH0tft)(0t0tft p间接法。它的特点是,在每一步迭代中都要满足
3、取极小的必要条件,而且要同时积分状态方程和协态方程,两种方程的积分都从 到 或从 到 。常用的间接法有边界迭代法和拟线性化法。H0t0tftft7.1 直接法直接法(一)梯度法。这是一种直接方法,应用比较广泛。它的特点是:先猜测任意一个控制函数 ,它可能并不满足 取极小的必要条件,然后用迭代算法根据 梯度减小的方向来改善 ,使它最后满足必要条件。)(tUH)(tUH 计算步骤如下:1 先猜测 中的一个控制向量 ,是迭代步数,初始时 。的决定要凭工程经验,猜得合理,计算收敛得就快。,0ftt)()(0tUtUK0K0UK2 在第 步,以估计值 和给定的初始条件 ,从 到 顺向积分状态方程,求出状
4、态向量 。KKU)(0tX0tft)(tXK 3 用 、和横截条件求得的终端值 ,从 到 反向积分协态方程,求出协态向量 。)(tK0tft)(ft)(tXK)(tUK4 表示在 、处取值。当这些量非最优值时,。0KgKKXKUKUH)(KKUHg)(计算哈密顿函数 对 的梯度向量KgUH5、是一个步长因子,它是待定的数。选择 使指标达到极小。这是一维寻优问题,有很多现成的优化方法可用。如分数法,0.618法,抛物线法,立方近似法等。(7-3)表明迭代是沿着梯度 的负方向进行的。KKKgKKKKgUU1修正控制向量(7-3)6、)()()(1KKKUJUJUJ计算是否满足下列指标(7-4)1
5、KK 是指定小量,若满足则停止计算,否则,令 ,转步骤2。另一停止计算的标准是Kg(7-5)例7-1考虑下面的一阶非线性状态方程uxx210)0(x(7-6)用梯度法寻找最优控制使下面的指标最小1022)(21dtuxJ(7-7)解因 自由,由横截条件得)1(x0)1(0哈密顿函数为uxuxH222)(21(7-8)协态方程为xxxH2(7-9)、选初始估计 。0)(0tu代入初始条件:,确定积分常数10)0(,0 xt101c、将 代入状态方程(7-6)可得0)(0tudtxdx2(7-11)ctx1积分上式可得(7-12)11010)()(0ttxtx代入(7-12)式即可得(7-13)3
6、将 代入协态方程(7-9),且由边界条件 从t=1倒向积分可得)(0tx0)1(00)1(121/)101(1 21)(020tt5 。这里选步长因子 。如此继续下去,直至指标函数随迭代变化很小为止。121/)101(1 21)()()(2001tuHtutu1K)()(00tuHuuH 4由 图7-1用梯度法寻找最优控制 图7-2 最优状态的求解01)(1tu最优值)(0tuut10和最优值)(tx)(0txxt 图7-1和图7-2表示了控制和状态的初始值和第一次迭代值,可以看到第一次迭代 就几乎收敛到最优值,与最优值还有差异,而且一般说来愈接近最优值收敛愈慢)(tx)(1tu梯度法应用得比
7、较多,它的优点是:(1)简单,编制程序容易;(2)计算稳定可靠。缺点是:(1)在接近最优解时,迭代收敛很慢,为改善 收敛性可用共轭梯度法和二阶变分法等;(2)不能区分局部极小和全局极小;(3)对控制变量受约束,终端状态受约束的情 况不能直接处理。对于这种有约束的情况可用约束梯度法或惩罚函数法加以处理。mibtuaiii,2,1)(约束梯度法可处理如下的不等式约束:(7-15)首先,对于任何控制 ,定义约束算子uCu iiiiiiiiiiiuibtubtuaatubtuatuCtu)()()()()()((7-16)显然 满足约束,即mitui2,1,)(满足约束,其中 ,再由 用无约束的梯度法
8、求解,在每一次迭代中得出 ,然后用 代替,再进行下一次迭代。u Tmuuu1u uCuuTmuuu1uCuu(7-17)惩罚函数法可处理如下形式的约束:qituXgi,2,10),((7-18)littXhffi,2,10),((7-19)这时,将性能指标 增广为)(uJ).(uJ其中,0,0iN0010)(iiiggg(7-21)fttqiliffiiiittXhNdtgtUXguJuJ01122),()(),()()((7-20)显然,当满足约束时,中后两项为零。当不满足约束时,后两项将使 增大,故称为惩罚函数。在迭代过程中,逐次增大 和 。显然当 和 很大时,所求得的 的无约束最优控制近
9、似于 的有约束最优控制。)(uJ)(uJiNiN)(uJ)(uJ(二)共轭梯度法 用共轭梯度法寻找最优控制时是沿着所谓共轭梯度向量的方向进行的。为了说明共轭梯度的意义,我们先从求函数极值问题的共轭梯度法开始,再推广到求泛函极值问题。1求函数极值的共轭梯度法其中,),(,)(2121nTTnaaaaxxxX为常数,为正定阵。CQ是 和 的内积。要求寻找 使 取极值。XQXX)(XFQXXQXXT),(7-23)设 是定义在 空间中的二次指标函数)(XFnRCXaQXXXFT),(21)((7-22)定义 则称 和 是 共轭的。(单位阵)时,共轭就变为通常的正交。XYQIQ 若 中两个向量 和 满
10、足nRXY0),(QYXQYXT(7-24)设向量 ,是两两 共轭的,以 为寻找方向,可得共轭梯度法的迭代寻优程序:KP,2,1,0KQKPKKKKPXX1(7-25)与梯度法不同处仅在于用共轭梯度 代替负梯度 。问题是如何产生共轭梯度方向 。,2,1,0,KPKKKXFg)(KP 值由 和 对 共轭的关系来确定,即QKPK1KP1KKKKPgP(7-26)0),(1KKQPP(7-27)令 ,即初始时共轭梯度与梯度方向相反、大小相等。以后的共轭梯度可如下递归产生:00Pg),(),(0111KKKKKKQPPgQPP),(),(111KKKKKQPPQPg将(7-26)代入(7-27),得K
11、称为共轭系数。),(),(111KKKKKQPPQPg故(7-28)用(7-28)式计算 是不方便的,因为要用到二阶导数阵 。由(7-22)和(7-23)知KQ 分别为 的第 个和第 个分量,右端表示由 的第 行第 列元素构成的矩阵。计算这个二阶导数阵非常困难。为此,有必要推导不用 来计算 的公式。jixx,XijQijQKnjixxXFQji,1)(2(7-29)性质1 若 是 空间中彼此 共轭的向量,则它们是线性独立的。110,nPPPnRQ在这个推导中要用到共轭梯度的下列性质:证明:证明:因为 正定,上式对每一个 成立,所以必须有 与假设矛盾,这说明 是线性独立的,它们构成了 空间中的一
12、组基向量。QjP1,2,1,0,0njCj110,nPPPnR上式左端各项对 取内积后有jQP100),(),(nijjjjiiQPPCQPPC(7-31)用反证法。若不线性独立,则必存在不全为零的常数 使 121,nCCC100niiiPC(7-30)其中,可这样来求:作内积1,1,0,niCi),(),),(10KKKKiniiKPQPCPPQCPQX(7-33)1210),(),(nKPQPPQXCKKKK,从而按照这个性质,函数 的极小点 可用这组基来表示,即)(XF*XX 10niiiPCX(7-32)性质2 式中,。(7-34)说明,在 处函数 的梯度 与前一步的寻找方向 必正交。
13、KXXKXFXF)()(1KXX)(XF1KgKP如果 ,则有0)(KKXFg0),(1KKgP(7-34)若不然,不妨先设 。再设 ,即 是最优步长。在 附近选一个 ,将 在 处展开,保留一阶项后,有 0),(1KKgP)(min0KKPXFKKKK0)(XFKKKKPXXX1)(KKKPXF证明:证明:KKKTKKKKKKKPXFPXFPXF)()()()(0100),)()(1010KKKKKTKKKPgPg(7-35)这与 为极小相矛盾。)(KKKPXF若设 则可取 ,同样得出矛盾,于是必有(7-34)成立。0),(1KKgPKK0性质3(7-36)说明,在 处的梯度 与以前各步的共轭
14、梯度寻找方向都正交。)(XFKXX Kg若 ,则必有0)(KKXFgKjPgjK10),(1(7-36)证明证明由(7-22)式所假定的二次函数 ,可得)(XFKKKQXaXFg)((7-38)101KjXXKjiiijK(7-37)得到111KKKKPXX重复使用设 为极小点,则*XX 0)(QXaXg(7-39)(7-38)减去(7-39)得)(XXQgKK(7-40)上式两边对 作内积,得1jP(7-42)),(),(),(1111jiKjiijjjKPQPPgPg1KjiiijKQXPQQXg101KjPQgKjiiij(7-41)=(7-37)代入(7-40),得由性质2知 再由 与
15、是 共轭的定义可知(7-42)右端第二项也为零,,0),(1jjPg1,kjiPi1jPQKjPgjK10),(1因此 (7-36)得证。但 是线性无关的,它们构成 中一组基,与所有基正交,而 中只有 个基,故 。这说明 处的梯度为零,即 为二次函数 的极小点。110,nPPPnRngnRn0ngnXnX)(XF如果取 ,则nK njPgjn10),(1(7-43)如果一个算法能在有限步内求出二次函数的极小点,就称这个算法具有二阶收敛性或有限步收敛性。由此可见,在 空间中,对二次函数 用(7-25)式所示的共轭梯度法寻优,迭代至多 步就可达到极小点。nR)(XFn性质4若 ,则0)(KKXFg
16、,2,10),(1KggKK(7-44)证:证:),(),(),(012111KKKKKKKKggPggPg0),(),(121KKKKKggPg(7-44)得证。由性质3和(7-26)式知 下面根据这四个性质来推出 的一个简单的计算公式。在(7-41)中令 ,可导出K1 Kj111,(),(KKKKKKQPgPgP),(),(),(1111KKKKKKKgPQPPgP),(11KKKKgPg),(),(111KKKKKgPgg),(),(),(),(),(),(12111112111KKKKKKKKKKKKKKKKKKKgPgggPgggPggPg再利用(7-26)式,可得 由性质4知 ,因
17、此得0),(1KKgg),(),(11KKKKKgPgP(7-45)用(7-46)计算 ,只用到 在 和 两处的梯度,因此非常方便。K)(XFKX1KX (7-46)对二次函数是精确的,对非二次函数,它只是一个近似公式 由性质3,就可得出21211),(),(KKKKKKKgggggg(7-46)将共轭梯度法求 的极小解的算式归纳如下:)(XF(4)递推逼近极值点解 用一维寻优决定。KKKKPXX1K212KKKgg00(2)算共轭系数 ,KKXFg)((1)算梯度001,gPPgPKKKK(3)算共轭梯度2、用共轭梯度法解最优控制问题 前面已说过,最优控制计算的直接法是用迭代方法逐步改善控制
18、量 ,使它最后满足哈密顿函数 取极小的必要条件。)(tuH 除了这些以外,其它在形式上与求函数极值的共轭梯度法一样。故梯度向量为KtutuKKuHuHtggK)()()()()((7-47)这里梯度向量 是时间的函数,向量时间函数的内积定义为)(tgKfttKKTKKKtgdttgtgtgtg02)()()()(),((7-48)共轭梯度法求最优控制步骤为(1)(2)(3)设已求出第K步估计的控制函数 可任选。)(),(0tutuK以 为初值,从 到 积分状态方程,得出状态轨迹 。)(0tX0tft)(tXK以 为终值,从 到 反向积分协态方程,求得协态轨迹 。)(ftft0t)(tK(4)(
19、5)(6)计算梯度向量kuuKuHg)(计算共轭系数时,。(7-49)212)()(tgtgKKK0K00时,。(7-50)1KKKKPgP0K00gP 计算共轭梯度(7)(8)停止计算。否则令 ,回到步骤2。1 KK当满足下面的不等式)()()(1KKKuJuJuJ(7-53)用一维寻优决定 ,即 KKKKKKPuJPuJ0min(7-52)KKKKPtutu)()(1(7-51)计算控制函数例7-2要求用共轭梯度法决定最优控制 ,使 最小。)(tuJ性能指标)1(2xJ(7-56)设系统状态方程为21)0(11xux(7-54)0)0(21)1(2212xuuxux(7-55)解ctxH)
20、(0222(常数)(7-59))1(211uxH(7-58)协态方程为212121)1(uuxuuH(7-57)哈密顿函数为0)1()1(11u(7-62)1)(2t(7-63)故协态方程化为11)1()1(22cxJ(7-61)0)1()1()1(111xJx(7-60)由横截条件(1)选 ,代入状态方程和协态方程(7-54)、(7-55)、(7-62)和(7-63),0K0)(0tu时的计算1,0221121xxx可求得1)(,21)(,21)(010201ttttxtx积分可得 02212100)()()(uxuHtgttx25121102010201梯度向量2500tgP共轭梯度 。(
21、2)时 ,用一维寻优来决定。将 代入状态方程(7-54)、(7-55)和协态方程(7-62)、(7-63),得1K)25()()(00001tPtutu01u)25()(01ttx2000102)25(21)25()25(1)(ttxttx21)252()(2011tttx)82581512138()()41526(2)(2342023012tttttttttx2002)(2497124921)1(xJ01297124900J积分得可求得 的最优值为097490)25(9749)(1ttu于是 )25(97491 1t由(7-62)式9711945119449)(211ttt1)(12t122
22、12111)()()(uxuHtg97229748t积分上式可得 210210220211194192)25()97229748()()(dttdtttgtg共轭系数2220111194901619418816)25(19419297229748tttPgP共轭梯度(3))194901619418816()25(9749)()(2211112ttPtutu2K)(2tu时,控制量为同以上步骤,将 代入状态方程和协态方程,求出2()u t)1)()()(),(22212221tttxtx和对 寻优,可得 ,于是1196194123)(2 ttu)()1(12fxJ由 所以这个例子只要两步迭代即可
23、得到最优解。一般说来,共轭梯度法比梯度法收敛快,但接近最优解后收敛性仍是较慢的。一个补救办法是重新启动,即找出几个共轭梯度方向 后,令 ,再用(7-50)重新迭代,寻找共轭梯度方向。110,nPPPnngP可以证明 ,即为最优控制。这只要证明)()(2tutu0)()(2221212uxuH即可。7.2 间接法(一)边界迭代法这个方法的特点是逐步改善对缺少的初始条件的估计,以满足规定的边界条件。它的原理如下。0),(UtUXH可解出 ,将它表示为 和 的函数,即UX利用哈密顿函数H取极小的方法),(tXUU(7-64)将所求得的 代入正则方程(7-1),消去正则方程中的 。再引入增广状态),(
24、tXUUUnRtY2)()()()(ttXtY(7-65)设(7-65)式有 个已知初始条件 ,个终端条件已知,设为n00)(XtXn 和 ,这是混合式的两点边值条件(参见例3-6),用边界迭代法也很易处理。),1()(qitxfi),1()(nqitfig一般是非线性向量函数。则正则方程(7-1)可写成),(tYgY(7-66)显然,是已知的,设)(ftfft)(Tnqqtttxtxt)()(),(),()(11(7-67)定义因 未知,用一个估计值 得到的解为)(0t)(0t)()(0ttf(7-69)设由 、出发积分正则方程(7-66),求得解 ,从中抽出 个分量构成 。显然 的值将随
25、而变,记成)(0tX)(0t)(tYn)(tZ)(ft)(0t)()(0ttf(7-68)因 估计得不一定准确,故 一般不等于给定值 .将(7-68)在 处展开为台劳级数,保留一次项,得)(0t)(ftf)()(00tt 其中,是 维矩阵,称为敏感矩阵或转移矩阵。Tnn)()()()()(0000tttttTf(7-70)式中,是 的第 行,第 列元素。(7-71)式右端表示由第 行第 列元素构成的矩阵。由(7-69)和(7-70)可得)()(0ttjfiijij)()(0000)()()(ttjfiTjjttZt(7-71))()()()()(1000ffTttttt(7-72)因 一般是非
26、线性函数,(7-72)式是一个近似式,为了求得正确的 ,要用迭代求解。)(0t 其中,是迭代次数,是松驰因子,可改善收敛性,收敛到最后时,将 取为1。在第 步,用 作为估值,积分正则方程,求得 ,K10K)(0tK)(fKt 令 是第 步的估值,则根据(7-72)可得到下面的迭代格式)(0tKK)()()()()(10001fKfKTKKttttt(7-73)为指定的小值,则停止计算。否则用 代替 ,再积分正则方程,重复进行。)(0tK)(01tK)()(fKftt若(7-74)计算步骤如下(1)(2)由 解出 ,代入状态和协态方程。0uH),(tXuu设已求出 的第 步估计值 和给定的合在一
27、起,从 积分正则方程,求出 抽出 个要求的分量的终值 ,若 ,停止计算,否则进行下一步。)()(fkftt)(0tK)(0tK)(0tXftt 到0)(),(ttXKKn)(fKt(3)(4)(5)按(7-73)计算 。)(01tK令 回到步骤2。1 KKT求敏感矩阵这种方法的缺点是:(1)(2)(3)第一次估计 很困难,)(0t终端值对 非常敏感时,与 相差很大,线性关系(7-70)不成立。)(0t)(ft)(ft敏感矩阵难于确定得很精确,对它求逆的运算也容易引入误差。例7-3 系统状态方程为性能指标为用边界迭代法寻找 ,使 最小。)(tu)(uJ1.00221)()(duuxuJ(7-77
28、)5)0(121xxx 5)0(414.04.1232212xuxxxx(7-76)(7-75)解因终端 ,自由,故)1.0(1x)1.0(2x0)1.0()1.0(21 设 的初始估计值为零,迭代结果见表7-1。可见在第7次迭代时,、已为零,满足了边界条件。)(0201tt(和)(1ft)(2ft表7-1 这个方法的特点是用迭代算法来改善对正则方程解的估计,使它逐步逼近正则方程的精确解。和前面一样,将正则方程写成。(二)拟线性化法 设已知 个初始条件 和 个终端条件 n00)(XtXnfft)(XY(7-79)),(tYgY(7-78)拟线性化法将非线性两点边值问题转化为线性两点边值问题,因
29、此变得容易求解。设在迭代的第 步获得近似解 ,将正则方程(7-78)在展开,保留一次项,可得到 步的近似解 ,有K)(tYK)(tYK1K1kY1,0)()(),(11KYYYgtYgYKKKKK(7-80)满足给定边界条件njXtXtYjjKj2,1)()(0001(7-81)nnnjttYjffjfKj2,2,1)()(1(7-82)(7-80)可写成下面的线性非齐次方程)(),()(11KKKKKKYYgtYgYYgY(7-83))()(11tYtAYKKKK或(7-84)是 的系统矩阵,nn22 KKYgtA)()(其中(7-85)可停止计算当满足nitYtYiKiKi2,2,1)()
30、(1(7-87)是 驱动函数向量。(7-84)是线性微分方程,由给定的 个边界条件可确定其通解的 个未知常数,故解 可完全被确定。1nn2n2)(1tYK)(),()(KKKKYYgtYgt(7-86)例7-4用拟线性化法求 ,使 最小。)(tuJ系统方程为性能指标为1022)(21dtuxJ(7-89)10)0(2xuxx(7-88)解u哈密顿函数为)()(21222uxuxH(7-90)0uuH(7-91)上式代入状态方程后得到10)0(2xxx(7-92)0)1()1(2XxxXH(7-93)xxxxY22或写成(7-94)上式与(7-78)对照可知xxxggtYg2),(221(7-9
31、5)根据(7-85)、(7-86)可得)(),()(KKKKYYgtYgtKKKKKKKKKKxxxxxx212122)(2KKKKKKxxgxggxgYgtA21212)()(2211(7-96)KKKxx2)(2(7-97)于是线性化后的正则方程(7-84)中的系数阵 和驱动项 都已确定,解这个非齐次时变微分方程,并用边界条件 和 以决定通解中的未定常数,就完全确定了 ,这就完成了一次迭代。当满足(7-87)式时,停止计算,求解结束。)(tAK)(tK10)0(1Kx0)1(1K)(1tYK7.3 小结小结 1 最优控制的计算方法可分为直接法和间接法两大类。直接法中我们列举了梯度法和共轭梯
32、度法。间接法中列举了边界迭代法和拟线性化法。2 直接法的特点是:在每步迭代中并不满足哈密顿函数 取极小的必要条件,只是在迭代终了才满足这个条件;另外积分状态方程时是从 ,而积分协态方程时是从 。由于状态和协态的稳定性是相反的,所以这种双向积分,可使最优化过程非常稳定。这可举一简单例子来说明。Hftt 到00ttf到例 7-5 uaxx(7-98)fttdtuJ0221(7-99)解:)(212uaxuH(7-100)于是协态方程为axH(7-101)设 ,则从 ,收敛,而 发散。若从 ,则 就变为收敛了。0aftt 到0)(tx)(t0ttf到)(t解这两个方程,得dueetxtxtttatt
33、a)()()(00)()(0(7-102))(00)()(ttaett(7-103)3 梯度法是利用梯度信息 来不断改善对控制函数 的估计,最后满足 的必要条件。这是一种简单又稳定的算法,几乎对所有的 的初始估计都有很好的收敛性。但在远离最优解时收敛速度快,在接近最优解时收敛得慢(原因在于 )。uH)(tu0uH)(tu0uH 共轭梯度法比梯度法稍微复杂些,但收敛速度也快些。同样,在接近最优解时,共轭梯度法收敛速度变慢。要加速接近最优解时的收敛速度可用二阶变分法,不过这种方法的计算复杂程度要增加很多。4 间接法的特点是:在每步迭代中都满足 取极小的必要条件;另外,它同时从一个方向(从 或从 )
34、积分状态和协态方程。Hftt 到00ttf到 由于状态和协态的稳定性相反,这就使得对边界条件的初始估计非常敏感。尤其当终端时刻远远大于系统的最小时间常数时,收敛性可能很差。5 边界迭代法是在每步迭代中不断改善对缺少的初始条件的估计去满足给定的终端条件(也可改善对缺少的终端条件的估计去满足给定的初始条件),这种方法对初始估计(如 )是非常敏感的,只有在能获得良好的初始估计时,才建议使用这种方法。)(0t 拟线性化法将非线性正则方程围绕上一步的估计解轨迹 线性化,递推解出 最后满足正则方程。这种方法对初始估计 可能不如边界迭代法对初始估计 那么敏感。并且求线性微分方程的解也比较容易。)(tYK),(1tYK)(0tY)(0t