一种工业加热炉系统的分数阶动态矩阵控制方法
【专利摘要】本发明公开了一种工业加热炉系统的分数阶动态矩阵控制方法。本发明首先采用Oustaloup近似方法将分数阶模型近似为整数阶高阶模型，基于近似高阶模型实施阶跃响应实验，采集阶跃响应数据，获取模型向量，然后将整数阶动态矩阵控制方法扩展到分数阶动态矩阵控制方法中，将分数阶微积分算子引入目标函数，进而基于阶跃响应模型和选取的目标函数设计了分数阶动态矩阵控制器。本发明运用于分数阶模型描述的实际过程对象，改善了整数阶DMC方法控制分数阶系统的不足之处，同时增加了调节控制器参数的自由度，获得了良好的控制性能，并能很好地满足实际工业过程的需要。
【专利说明】
一种工业加热炉系统的分数阶动态矩阵控制方法
技术领域
[0001] 本发明属于自动化技术领域，涉及一种工业加热炉系统的分数阶动态矩阵控制 (FDMC)方法。
【背景技术】
[0002] 在实际工业控制过程中，随着对产品的控制精度和安全操作的要求越来越高，许 多复杂的对象是整数阶微分方程无法精确描述的，用分数阶微分方程能更精确地描述对象 特征和评估产品性能。动态矩阵控制(DMC)方法作为先进控制方法的一种，具有对模型要求 低，鲁棒性强，处理延时的方法简单易行等特点，在实际过程控制中获得大量成功的应用。 针对一类单输入单输出分数阶系统，传统PID控制方法和整数阶DMC方法对这类对象的控制 效果并不是很好，这就需要研究具备良好控制性能的控制器来控制这类用分数阶模型描述 的实际被控对象。如果将整数阶DMC方法扩展到分数阶DMC方法中，那将能有效的弥补整数 阶DMC方法在控制分数阶系统的不足，并能获得更好的控制效果，同时也能促进DMC在分数 阶系统中的运用。
【发明内容】
[0003] 本发明的目的是针对分数阶模型描述的加热炉温度对象，提供一种工业加热炉系 统的FDMC方法，以维持分数阶系统的稳定性并保障良好的控制性能。该方法首先采用 Oustaloup近似方法将分数阶模型近似为整数阶高阶模型，基于近似高阶模型实施阶跃响 应实验，采集阶跃响应数据，获取模型向量，然后将整数阶动态矩阵控制方法扩展到分数阶 动态矩阵控制中，将分数阶微积分算子引入目标函数，进而基于阶跃响应模型和选取的目 标函数设计了分数阶动态矩阵控制器。
[0004] 该方法可以很好地运用于分数阶模型描述的实际过程对象，改善了整数阶DMC方 法控制分数阶系统的不足之处，同时增加了控制器参数调节的自由度，获得了良好的控制 性能，并能很好地满足实际工业过程的需要。
[0005] 本发明的技术方案是通过数据采集、模型建立、预测机理、优化等手段，确立了一 种工业加热炉系统的FDMC方法，该方法可有效提高系统的控制性能。
[0006] 本发明方法的步骤包括：
[0007] 步骤1、建立加热炉中温度对象的分数阶阶跃响应模型，具体方法是：
[0008] 1.1采集实际过程对象的实时输入输出数据，利用该数据建立被控对象在t时刻的 分数阶微分方程模型，形式如下：
[0009] (〇 + qyai)(f) + C〇少⑴二.M.(.'f)
[0010] 其中，ai，a2为微分阶次，c〇，Cl， C2为相应的系数，y(t)，u⑴分别为过程的输出和输 入。
[0011] 1.2根据分数阶微积分定义，对步骤1.1中的模型进行拉氏变换，得到被控对象的 传递函数形式如下：
[0013]其中s为拉普拉斯变换算子。
[0014] 1.3由Oustaloup近似方法得到微分算子8°的近似表达形式如下：
[0016] 其中，a为分数阶微分阶次，0〈a〈 1，N为选定的近似阶次，[=<， Wn = wbwu~" ^ _<X VN > wn = WbW(^l+a }/N , wu =」wh/wb，Wb 和 Wh 分别为选定的拟合频率 的下限和上限。
[0017] 1.4根据步骤1.3中的方法，将步骤1.2中的分数阶系统近似为整数阶高阶系统，给 所得高阶模型一个阶跃输入信号，记录高阶模型的阶跃响应曲线。
[0018] 1.5将步骤1.4得到的阶跃响应曲线进行滤波处理，然后拟合成一条光滑曲线，记 录光滑曲线上的每个采样时刻对应的阶跃响应数据，第一个采样时刻为T s，相邻两个采样 时刻间隔的时间为Ts，采样时刻顺序为Ts、2T S、3TS……；高阶模型的阶跃响应将在某一个时 亥lJtN = NTs后趋于平稳，当at(t>N)与aN的误差和测量误差有相同的数量级时，即可认为 似等于阶跃响应的稳态值。建立高阶模型的模型向量a:
[0019] a=[ai，a2,…，aN]T
[0020] 其中T为矩阵的转置符号，N为建模时域。
[0021 ]步骤2、设计被控对象的分数阶动态矩阵控制器，具体方法如下：
[0022] 2.1利用步骤1获得的模型向量a建立被控对象的动态矩阵，其形式如下： ctx 0 * " 0 Clr, U] ? ? 0
[0023] A-； , . * * ? * * * ? nT ap . * aP
[0024] 其中，A是被控对象的PXM阶动态矩阵，ai是阶跃响应的数据，P、M分别为动态矩阵 控制算法的优化时域和控制时域，M〈P〈N。
[0025] 2.2求取被控对象当前k时刻的模型预测初始响应值yM(k)
[0026] 首先，在k_l时刻加入控制增量A u(k_l)后得到模型预测值yp(k_l):
[0027] yp(k_l) =yM(k_l )+Ao A u(k_l)
[0028] 其中，
[0029] yp(k_l) = [yi(k | k_l)，yi(k+l | k_l)，…，yi(k+N_l | k_l) ]T
[0030] yM(k-l) = [yo(k | k_l)，yo(k+l | k_l)，…，yo(k+N_l | k_l) ]T
[0031 ] Ao= [ai，a2,…，aN]T
[0032] yi(k | k_l)，yi(k+l | k_l)，…，yi(k+N_l | k_l)分别表示被控对象在k_l时刻对k，k+ 1，…，k+N-1 时刻的模型预测值，yo(k|k-l)，yo(k+l |k_l)，…，yi，o(k+N-l |k_l)表示k-1 时刻 对k，k+l，…，k+N-1时刻的初始预测值，Ao为阶跃响应数据建立的矩阵，Au(k-l)为k-1时刻 的输入控制量；
[0033] 然后，可以得到k时刻被控对象的模型预测误差值e(k):
[0034] e(k) = y(k)-yi(k | k~l)
[0035] 其中y(k)表示k时刻测得的被控对象的实际输出值；
[0036] 进一步得到k时刻修正后的模型输出值yc^r(k):
[0037] ycor(k) =yM(k_l)+h*e(k)
[0038] 其中，
[0039] yc〇r(k) = [yc0r(k | k)，yc0r(k+l | k)，…，yc0r(k+N_l | k) ]T
[0040] h=[l，a，…，a]T ] yc〇r(k | k)，yc〇r(k+l | k)，…，yc〇r(k+N_l | k)分别表示被控对象在k时刻模型的修正 值，h为误差补偿的权矩阵，ct为误差校正系数；
[0042] 最后得到k时刻的模型预测的初始响应值yM(k):
[0043] yM(k) = Sy cor (k)
[0044] 其中，S为NXN阶的状态转移矩阵， 0 1 〇 0 G 1 0 …
[0045] S =：： ： ' … 0 0 0 …1 0 0 …0 1
[0046] 2.3计算被控对象在M个连续的控制增量A U(k)，A U(k+1)，…，A u(k+M-l)下的预 测输出值y?，具体方法是：
[0047] ypM(k) =yp〇(k)+A A UM(k)
[0048] 其中，
[0049] ypM(k) = [yM(k+l |k)，yM(k+2|k)，…，yM(k+P|k)]T
[0050] yp〇(k) = [yo(k+l |k) ,yo(k+2|k) ,**? ,yo(k+P|k)]T
[0051] AuM(k) = [ Au(k)，Au(k+1)，."，Au(k+M_1)]T
[0052] yp〇(k)是yM(k)的前P项，yM(k+l | k)，yM(k+2 | k)，…，yM(k+P | k)为k时刻对k+1，k+ 2,…，k+P时刻的模型预测输出值。
[0053] 2.4选取被控对象的参考轨迹和动态矩阵控制方法的目标函数JFDMC，其形式如下：
[0054] yr(k+i)=入1yp(k) + ( 1_入Odk) J郷ft：. = '1 A;/、[)',.⑴-⑴]一 + ': 乂7::《卜 1 疒
[0055] .prs tMTs , ^ =...!'[yr(t)-y^)Ydt+ f o' ^-Au(i-\) dt
[0056] 其中，y 1, y 2为任意实数，7/p/0)表示函数f(t)在[ti,t2]上的y次积分，D为微 分符号。
[0057] 依据Griinwald-Letnikov分数阶微积分定义，对上述目标函数在采样时间Ts进行 离散化，得到：
[0058] 三[y;.⑷-v.w (々)]' A()V 7:)[y;.⑷->，? ^
[0059] 其中，
[0060] Yr(k) = [yr(k+l),yr(k+2),??? ,yr(k+P) ]T
[0061 ]八(Ye, I；) = ，Wp-2，…,Wi，w0)
[0062] W =0)(Ys) -（0(Ys> 9 q ?-(卢一1)
[0064]在上式中进一步引入误差加权系数Q = diag(qi，q2，???，qp)和控制加权系数R = diag(ri，r2，'"，rp)，所得目标函数为
[0065] JFDMC=[Yr(k)-ypM(k)]TA ( y i,Ts)Q[Yr(k)-ypM(k)]+AUTA ( y2,Ts)RAU
[0066] 2.5依据步骤2.4中的目标函数求解
，得到控制量，形式如下：
[0067] AUM(k) = (ATA( yi,Ts)QA+A( y2,Ts)R)_1ATA( yi,Ts)Q(Yr(k)-yp〇(k))
[0068] A u(k) = [ 1,0，…，0] A UM(k)
[0069] u(k) =u(k-l)+ A u(k)
[0070] 2.6在k+1时刻，1 = 1，2,3,…，依照2.1到2.5中的步骤依次循环求解分数阶动态矩 阵控制器的控制量u(k+l)，再将其作用于被控对象。
[0071 ]本发明提出了一种工业加热炉系统的分数阶动态矩阵控制方法，该方法建立了被 控对象的分数阶阶跃响应模型，将微分算子引入控制器增加了控制器参数调节的自由度， 有效弥补了整数阶动态矩阵控制对分数阶系统控制的不足之处，提高了系统的控制性能， 同时促进了动态矩阵控制方法在分数阶系统中的应用。
【具体实施方式】
[0072] 以实际过程中加热炉的温度过程控制为例：
[0073] 由加热炉的实时温度数据得到分数阶模型，温度控制系统的调节手段是控制阀门 开度。
[0074] 步骤1、建立加热炉中温度对象的分数阶阶跃响应模型，具体方法是：
[0075] 1. 1采集加热炉温度对象的实时输入输出数据，利用该温度数据建立加热炉温度 对象在t时刻的分数阶微分方程模型，形式如下：
[0076] C2^(K2)(〇 + 6^(〇1)(〇 + <?〇>'(〇 = li(t)
[0077] 其中，(^，句为微分阶次，CQ，C1，C2为相应的系数，y(t)，u(t)分别为加热炉温度对象 的温度输出和阀门开度。
[0078] 1.2根据分数阶微积分定义，对步骤1.1中的模型进行拉氏变换，得到温度对象的 传递函数形式如下：
[0080]其中S为拉普拉斯变换算子。
[0081 ] 1.3由Oustaloup近似方法得到微分算子8°的近似表达形式如下：
[0083] 其中，a为分数阶微分阶次，0〈a〈 1，N为选定的近似阶次，， < =WAWp-卜a)/A/，二 "―1+c〇W，=」Wh/wb，w b和wh分别为选定的拟合频率 的下限和上限。
[0084] 1.4根据步骤1.3中的方法，将步骤1.2中的分数阶系统近似为整数阶高阶系统模 型，给所得高阶模型一个阶跃输入信号，记录高阶模型的阶跃响应曲线。
[0085] 1.5将步骤1.4得到的阶跃响应曲线进行滤波处理，然后拟合成一条光滑曲线，记 录光滑曲线上的每个采样时刻对应的阶跃响应数据，第一个采样时刻为T s，相邻两个采样 时刻间隔的时间为Ts，采样时刻顺序为Ts、2T S、3TS……；高阶模型的阶跃响应将在某一个时 亥lJtN = NTs后趋于平稳，当at(t>N)与aN的误差和测量误差有相同的数量级时，即可认为 似等于阶跃响应的稳态值。建立高阶模型的模型向量a:
[0086] a=[ai，a2,…，aN]T
[0087]其中T为矩阵的转置符号，N为建模时域。
[0088] 步骤2、设计加热炉温度对象的分数阶动态矩阵控制器，具体方法如下：
[0089] 2.1利用步骤1获得的模型向量a建立加热炉温度对象的动态矩阵，其形式如下： 0J 0 ? * 0 a) 〇, ? ? 0
[0090] . jXp
[0091] 其中，A是温度对象的PXM阶动态矩阵，ai是阶跃响应的数据，P、M分别为动态矩阵 控制算法的优化时域和控制时域，M〈P〈N。
[0092] 2.2求取加热炉温度对象当前k时刻的模型预测初始响应值yM(k)
[0093] 首先，在k-1时刻加入阀门开度增量Au(k-l)后得到模型预测值yP(k-l):
[0094] yp(k_l) =yM(k_l )+Ao A u(k_l)
[0095] 其中，
[0096] yp(k_l) = [yi(k|k_l)，yi(k+l|k_l)，???，yi(k+N_l|k_l)]T
[0097] yM(k-l) = [yo(k | k_l)，yo(k+l | k_l)，…，yo(k+N_l | k_l) ]T
[0098] Ao= [ai,a2, ?" )&n]t
[0099] yi(k | k_l) ,yi(k+l | k_l),…，yi(k+N_l | k_l)分别表示温度对象在k_l时刻对k,k+ 1，???，k+N_l时刻的模型预测值，yo(k|k_l) ,yo(k+l |k_l),…，yi,o(k+N_l |k-l)表示k_l时刻 对k，k+1，…，k+N-1时刻的初始预测值，Ao为阶跃响应数据建立的矩阵，A u(k_l)为k-1时刻 的阀门开度增量；
[0100] 然后，可以得到k时刻加热炉温度对象的模型预测误差值e(k):
[0101] e(k) = y(k)-yi(k | k~l)
[0102] 其中y(k)表示k时刻测得的加热炉温度对象的实际输出值；
[0103] 进一步得到k时刻修正后的模型输出值ycOT(k):
[0104] yc〇r(k) =yM(k_l )+h*e(k)
[0105] 其中，
[0106] ycor(k) = [ycor(k | k)，ycor(k+l | k)，…，ycor(k+N-l | k) ]T
[0107] h=[l,a,---,a]T
[0108] | k) jc^k+l | k),…，ycxxrU+N-l | k)分别表示加热炉温度对象在k时刻模型 的修正值，h为误差补偿的权矩阵，a为误差校正系数；
[0109] 最后得到k时刻的模型预测的初始响应值yM(k):
[0110] yM(k) = Sy cor (k)
[0111] 其中，S为NXN阶的状态转移矩阵， 0 1 〇 ...... 0 0 1 0
[0112] S.二：? ： 0 0 0 1 _00???0 1 _
[0113] 2.3计算温度过程在M个连续的阀门开度增量A U(k)，A U(k+1)，…，A U(k+M-1)下 的预测输出值y?，具体方法是：
[0114] ypM(k)=yp〇(k)+AAuM(k)
[0115] 其中，
[0116] ypM(k) = [yM(k+l |k)，yM(k+2|k)，???，yM(k+P|k)]T
[0117] yp〇(k) = [yo(k+l |k) ,yo(k+2 | k), ??? ,yo(k+P | k)]T
[0118] AuM(k) = [Au(k), Au(k+1)，."，A u(k+M_l) ]T
[0119] yp〇(k)是yM(k)的前P项，yM(k+l | k)，yM(k+2 | k)，…，yM(k+P | k)为k时刻对k+1，k+ 2，…，k+P时刻的模型预测输出值。
[0120] 2.4选取温度对象的参考轨迹和动态矩阵控制方法的目标函数JFDMC，其形式如下： [0121 ] yr(k+i) =A1yp(k) + ( l-A^cCk) l,>u< =71 ^ [yr (n - y.u (of+v： ^ -if
[0122] c：mts T 7 =I ' 〇 I D] -Au(t~\Ydt 丨Ts: .打s
[0123]其中，丫 1，y 2为任意实数，7/^/⑴表示函数f(t)在[ti，t2]上的y次积分，D为微 分符号。
[0124] 依据Griinwald-Letnikov分数阶微积分定义，在采样时间Ts下对上述目标函数进 行离散化，得到：
[0125] ^//)>? = [}/- (k) ~ yp.M (^)]' A(r i ^ Ts)[}/- (k) ~ y；>M (/i：)]+Aiim (k)! A(y 2^ )Aum (k)
[0126] 其中，
[0127] Yr(k) = [yr(k+l)，yr(k+2)，…，yr(k+P) ]T
[0128] 八(YE,l) = 命g(W"丨,W" Wj, w(,)
[0129] W -03(Ts)
[0131 ]在上式中进一步引入误差加权系数Q = diag(qi ,q2, ??? ,qp)和控制加权系数R = diag(ri,r2,'"，rp)，所得目标函数为
[0132] JFDMC=[Yr(k)-ypM(k)]TA ( y i,Ts)Q[Yr(k)-ypM(k)]+AUTA ( y2,Ts)RAU
[0133] 2.5依据步骤2.4中的目标函数求解
，得到阀门开度增量，形式如下：
[0134] AUM(k) = (ATA( yi,Ts)QA+A( y2,Ts)R)_1ATA( yi,Ts)Q(Yr(k)-yp〇(k))
[0135] A u(k) = [ 1,0，…，0] A UM(k)
[0136] u(k) =u(k_l)+A u(k)
[0137] 2.6在k+1时刻，1 = 1，2,3,…，依照2.1到2.5中的步骤依次循环求解分数阶动态矩 阵控制方法的阀门开度u(k+l)，再将其作用于加热炉温度对象。
【主权项】
1. 一种工业加热炉系统的分数阶动态矩阵控制方法，其特征在于该方法包括以下步 骤； 步骤1、建立加热炉中温度对象的分数阶阶跃响应模型，具体是： 1.1采集实际过程对象的实时输入输出数据，利用该数据建立被控对象在t时刻的分数 阶微分方程模型，形式加下·其中，Ci1，α2为微分阶次，CO，C1，C2为相应的系数，y (t)，U(t)分别为过程的输出和输入； 1.2根据分数阶微积分定义，对步骤1.1中的模型进行拉氏变换，得到被控对象的传递 函数形式如下：其中s为拉普拉斯变换算子； 1.3由Oustaloup近似方法得到微分算子8°的近似表达形式如下：其中，α为分数阶微分阶次，〇〈 a〈 I，N为选定的近似阶次，和Wh分别为选定 的拟合频率的下限和上限； 1.4根据步骤1.3中的方法，将步骤1.2中的分数阶系统近似为整数阶高阶系统，给所得 高阶模型一个阶跃输入信号，记录高阶模型的阶跃响应曲线； 1.5将步骤1.4得到的阶跃响应曲线进行滤波处理，然后拟合成一条光滑曲线，记录光 滑曲线上的每个采样时刻对应的阶跃响应数据，第一个采样时刻为Ts，相邻两个采样时刻 间隔的时间为T s，采样时刻顺序为Ts、2TS、3TS……；高阶模型的阶跃响应将在某一个时刻t N = NTs后趋于平稳，当at(t>N)与aN的误差和测量误差有相同的数量级时，即可认为a N近似等 于阶跃响应的稳态值;建立高阶模型的模型向量a: a= [ai，a2，…，aN]T 其中T为矩阵的转置符号，N为建模时域； 步骤2、设计被控对象的分数阶动态矩阵控制器，具体如下： 2.1利用步骤1获得的模型向量a建立被控对象的动态矩阵，其形式如下：其中，A是被控对象的PXM阶动态矩阵，ai是阶跃响应的数据，P、M分别为动态矩阵控制 算法的优化时域和控制时域，M〈P〈N; 2.2求取被控对象当前k时刻的模型预测初始响应值yM(k) 首先，在k-1时刻加入控制增量Au(k-l)后得到模型预测值yp(k-l): yp(k-l) =yM(k_l)+A〇A u(k_l) 其中， yp(k-l) = [yi(k|k_l)，yi(k+l |k_l)，···，yi(k+N_l |k_l)]T yM(k-l) = [y〇(k|k-l)，y〇(k+l |k-l)，···，y〇(k+N-l |k-l)]T Ao= [ai，a2，···，aN]T yi(k|k-l)，yi(k+l |k-l)，···，yi(k+N-l |k-l)分别表示被控对象在k-1时刻对k，k+l，…， k+N-1 时刻的模型预测值，y〇(k|k-l)，y〇(k+l |k-l)，…，yi,〇(k+N-l |k-l)表示k-1 时刻对k，k+ I，…，k+N-1时刻的初始预测值，Ao为阶跃响应数据建立的矩阵，Au(k-l)为k-1时刻的输入 控制量； 然后，得到k时刻被控对象的模型预测误差值e(k): e(k)=y(k)-yi(k|k-l) 其中y (k)表示k时刻测得的被控对象的实际输出值； 进一步得到k时刻修正后的模型输出值y?r(k): yc〇r(k) =yM(k_l)+h*e(k) 其中， y〇〇r(k) - [yc〇r(k I k)，yc〇r(k+ l|k)，^"，yc〇r(k+N_11 k)] h=[l,a,---,a]T yc〇：r(k |k)，yCC)r(k+l I k), ···，yc〇r(k+N-l I k)分别表示被控对象在k时刻模型的修正值，h 为误差补偿的权矩阵，a为误差校正系数； 最后得到k时刻的模型预测的初始响应值yM(k): yM(k) =Syc〇r(k) 其中，S为N X N阶的状态转秤$ _2.3计算被控对象在M个连续的控制增量Δ u(k)，Δ u(k+l)，···，Δ u(k+M-l)下的预测输 出值ypM，具体是： ypM(k) =yp〇(k)+A Δ UM(k) 其中， ypM(k) = [yM(k+l I k)，yM(k+2 I k)，…，yM(k+P I k) ]τ yp〇(k) = [yo(k+l | k) ,yo(k+2 | k), ··· ,y〇(k+P | k) ]T AuM(k) = [Au(k), Au(k+1)，···，Δ u(k+M_l) ]T yp〇(k)是yM(k)的前P项，yM(k+l I k)，yM(k+2 I k)，…，yM(k+P I k)为k时刻对k+1，k+2，…，k+ P时刻的模型预测输出值； 2.4选取被控对象的参考轨迹和动态矩阵控制方法的目标函数Jrac，其形式如下：其中，γ 1，γ 2为任意实娄 €示函数f (t)在[ti，t2]上的γ次积分，D为微分符 号； 依据Griinwald-Letnikov分数阶微积分定义，对上述目标函数在采样时间Ts进行离散 化，得到：在上式中进一步引入误差加权系数0 = (^8(91，92，'"仰)和控制加权系数1? = (^8(1'1， r2，'"，rP)，所得目标函数为 Jfdmc= [Yr(k)-ypM(k)]T Λ ( γ i ,Ts)Q[Yr(k)-ypM(k) ]+Δ UtA ( y2,Ts)RAU 2.5依据步骤2.4中的目标函数求彳得到控制量，形式如下： AuM(k) = (ATA ( yi,Ts)QA+A ( y2,Ts)R)_1ATA ( γ i,Ts)Q(Yr(k)-yp〇(k)) A u(k) = [1,0,…，0] Δ UM(k) u(k)=u(k_l)+Au(k) 2.6在k+1时刻，1 = 1,2,3，···，依照2.1到2.5中的步骤依次循环求解分数阶动态矩阵控 制器的控制量u(k+l)，再将其作用于被控对象。
【文档编号】G05B13/04GK105892296SQ201610307224
【公开日】2016年8月24日
【申请日】2016年5月11日
【发明人】汪大卫, 张日东, 张俊锋
【申请人】杭州电子科技大学