The Closed-loop Supply Chain Network Equilibrium with Products Lifetime in Multi-period Planning Model
-
摘要: 在考虑产品寿命次数的基础上，研究了由多个制造商和多个消费市场组成的闭环供应链网络的多期动态均衡问题．采用变分不等式和互补理论获得了制造商和消费市场的最优性条件，进而建立了整个闭环供应链网络的均衡模型．采用定步长的修正投影收缩算法对模型进行求解，通过算例将不同寿命次数结合有关参数对网络均衡状态的影响进行了对比分析，给出了政府有关部门和供应链网络中制造商作决策时需关注的关键问题，本文所得结论为开展相关领域的研究提供一定的借鉴．Abstract: An equilibrium problem of multi-period dynamic closed-loop supply chain network,which comprises manufacturers and demand markets,based on the product's lifetime is studied. Using variational inequality and complementary theory,the optimal conditions for manufacturers and demand markets are formulated,and thus the integrated closed-loop supply chain network equilibrium model is established. The equilibrium is obtained by employing modified project contraction algorithm with fixed step. Using numerical examples,the impacts of different product lifetime and related parameters on the supply chain network equilibrium results are compared and analyzed,and the solution to a critical problem that is encountered while making decisions in the relevant departments of government and manufacturers in a supply chain network is provided. The results of this study can help promote research in other related fields.
-
0. 引言
自抗扰控制(ADRC)将系统的外部干扰与内部动态不确定性视为总扰动，运用扩张状态观测器(extented statae observer，ESO)对总干扰进行实时观测，并通过反馈控制器进行补偿[1]。低阶的线性ADRC有着结构简单、鲁棒性能强、参数整定简单等优势[2]，在学界和工业界得到了广泛的关注。为实现低阶ADRC对时滞系统[3]和高阶系统[4]的控制，学界常采用输出预估、输入预估和输入补偿[5]等方法改进低阶ADRC。然而在实际工业过程中，被控对象往往为系统中各回路之间相互耦合的多变量系统[6]，当其中含有高阶大惯性环节时，低阶ADRC控制这类高阶多变量系统将迎来很大的挑战。
针对高阶大惯性对象，有研究者提出了一种利用高阶系统模型信息进行补偿的自抗扰控制策略[7]，实现了低阶ADRC对高阶大惯性系统的良好控制。而对于多变量系统的耦合问题，学者们提出了分散控制、集中控制和解耦控制三种策略[8]。分散控制将多变量系统中各回路之间的耦合视为总干扰，将多变量系统视为独立的单变量系统进行控制，有着易于设计、调整和维护的优点[9]，但是其牺牲了一定的系统输出性能，不适用于耦合强的多变量系统。集中控制策略则是以系统复杂度和计算量为代价，通过相应的算法换取控制系统的高性能。而解耦控制则是运用解耦器将被控对象各回路之间的耦合作用削弱或消除，使得被控对象解耦为对角占优或对角的形式，再对各自独立回路的单变量对象进行控制器设计以达到满意的控制效果[10]。目前已有的解耦方法包括理想解耦[11]、简化解耦[12]、正规解耦[10, 13]、前馈补偿解耦[14]和逆解耦[15-16]等。分散式ADRC将回路之间的耦合视为扰动并为每条回路单独设计ESO进行观测。有研究者证明，这种方法对多变量系统具有近似等效逆解耦的效果，并且相较于普通的解耦控制策略有着更好的鲁棒性[17]，同时解决了集中控制策略计算量和系统复杂度大的问题，更有利于工程应用。分散式ADRC在一定程度上解决了回路之间的耦合问题，但是其采用的ADRC结构通常是低阶的，对于高阶多变量系统其控制效果往往不佳。
针对此问题，本文将分散式ADRC与低阶补偿ADRC相结合，提出了一种分散式补偿ADRC策略，同时给出了其系统结构与参数整定公式。为探究其是否改善了低阶分散式控制器对于高阶多变量系统的控制效果，本文根据分散式补偿ADRC的系统结构，推导出了其闭环传递函数。在此基础上，本文利用多变量系统中的逆奈奎斯特阵列设计方法，总结出了一种多变量控制系统的稳定区域定量分析方法，并将其运用在对分散式补偿ADRC的分析当中。最后，利用系统仿真实验和蒙特卡洛实验对分散式补偿ADRC的控制效果和鲁棒性能进行了检验。
1. 算法介绍
1.1 自抗扰控制原理
考虑一个1阶系统：
˙y=g(t,y,d)+bu (1) 其中，u、y、d、b分别表示系统的输入、输出、外部扰动和输入增益，g(t，y，d)表示系统中时变、扰动、动态不确定性等的综合。实际中输入增益的值可能是未知的，可定义b0为输入增益的估计值，同时定义f=g+(b-b0)u为系统总扰动，包括系统的外部干扰和内部不确定性，则式(1)的模型可以转化为
˙y=f+b0u (2) 定义1阶系统的状态向量为x =[y]，其只包含一个分量。为精准估计系统的总扰动f，将其扩张为新的系统状态，令x2=f，所以系统的状态向量被扩张为x =[x1 x2]T=[y f]T。式(2)中的系统可描述为状态空间的形式：
\left\{\begin{array}{l} \dot{\boldsymbol{x}}=\boldsymbol{A} \boldsymbol{x}+\boldsymbol{B} u+\boldsymbol{E} f \\ y=\boldsymbol{C} \boldsymbol{x} \end{array}\right. (3) 其中，
\boldsymbol{A}=\left[\begin{array}{ll} 0 & 1 \\ 0 & 0 \end{array}\right], \boldsymbol{B}=\left[\begin{array}{l} b_0 \\ 0 \end{array}\right], \boldsymbol{E}=\left[\begin{array}{l} 0 \\ 1 \end{array}\right], \boldsymbol{C}=\left[\begin{array}{ll} 1 & 0 \end{array}\right] (4) 针对式(3)中的系统，1阶ESO的设计为
\dot{\boldsymbol{z}}=\boldsymbol{A} \boldsymbol{z}+\boldsymbol{B} u+\boldsymbol{L}\left(x_1-z_1\right) (5) 其中，z =[z1 z2]T为状态向量的观测值，L =[β1 β2]T为观测器的增益向量。当L中的参数β1和β2整定合适时，z中的观测值z1和z2分别跟踪系统输出y和总扰动f。
对系统进行补偿的状态反馈控制律(state feedback control law，SFCL)设计为
u=\boldsymbol{K}(\boldsymbol{r}-\boldsymbol{z}) / b_0 (6) 其中， \boldsymbol{r}=\left[\begin{array}{ll} r & \dot{r} \end{array}\right]^{\mathrm{T}}表示系统的输入参考信号，K =[kp 1]为控制器的增益向量，kp为控制器的反馈增益。因工业过程对象通常为常值调节系统，所以系统输入r可认为恒值，故 \dot{r}=0。
综上，1阶ADRC的结构如图 1所示。基于文[18]提出的带宽参数化方法，观测器增益设置为β1=2ωo和β2=ωo2，因此A-LC的所有特征值都位于-ωo处。同时，控制器增益确保矩阵 \widetilde{A}的所有特征值都位于-ωc处，其中 \widetilde{A}定义为
\left[\begin{array}{ll} \widetilde{\boldsymbol{A}} & 0 \\ 0 & 0 \end{array}\right]=\boldsymbol{A}-\boldsymbol{B} \boldsymbol{K} (7) 可以得到kp=ωc。因此，需整定的ADRC参数简化为输入增益的估计值b0、控制器带宽ωc和观测器带宽ωo。
1.2 补偿自抗扰控制器
在热工领域等实际应用中，被控对象常被识别成高阶大惯性系统：
G(s)=\frac{K}{(T s+1)^n} (8) 其中，K、 T和n分别为高阶对象的增益、时间常数和阶数。在实际工程中，应用最广泛的ADRC控制器通常是低阶的，而由于低阶ESO的阶次与高阶对象的阶次不匹配，所以其无法观测到被控对象所有阶次的状态信息，最终导致低阶ADRC对于高阶对象的控制效果不佳。为解决这一问题，针对n阶对象设计m阶补偿ADRC控制器。为了尽可能地观测到系统更多的状态信息，只需把观测到的系统补偿至与于ADRC的阶次一致即可，故补偿环节设置为
G_{\mathrm{cp}}(s)=\frac{1}{(T s+1)^{n-m}} (9) 虽然补偿后丢失了被控对象的部分高阶状态信息，但是由于观测到的系统与ESO的阶次相匹配了，ESO的观测精度仍然获得了提高，所以补偿控制器能够得到更好的控制效果。加入补偿环节Gcp后系统的结构如图 2所示，将其称为补偿ADRC。
补偿ADRC的参数整定方法为[7]：
1) 利用机理建模或系统辨识方法将被控过程近似为G(s)=K/(Ts+1)n的形式；
2) 选择所需的补偿ADRC控制器阶次m，设计补偿环节Gcp(s)=1/(Ts+1)n-m；
3) 计算各控制器参数：
\left\{\begin{array}{l} b_0=\frac{K}{T^m} \\ \omega_{\mathrm{c}}=\frac{1}{T} \\ \omega_{\mathrm{o}}=t \omega_{\mathrm{c}}, t \in[10, 100] \end{array}\right. (10) 1.3 分散式补偿自抗扰控制器
1.3.1 问题描述
在多变量系统中，也存在着传递函数矩阵元素为高阶大惯性对象的情况。设Y =[Y1，Y2，…，Yn]T为系统输出向量；R =[R1，R2，…，Rn]T为系统输入向量；G =[Gij]为传递函数矩阵，Gij表示输入Rj到输出Yi的传递函数。n×n阶多变量系统的输入输出与传递函数矩阵关系为
\left[\begin{array}{c} Y_1 \\ Y_2 \\ \vdots \\ Y_n \end{array}\right]=\left[\begin{array}{cccc} G_{11} & G_{12} & \cdots & G_{1 n} \\ G_{21} & G_{22} & \cdots & G_{2 n} \\ \vdots & \vdots & \ddots & \vdots \\ G_{n 1} & G_{n 2} & \cdots & G_{n n} \end{array}\right]\left[\begin{array}{c} R_1 \\ R_2 \\ \vdots \\ R_n \end{array}\right] (11) 其中，主对角传递函数Gij(i=j)均辨识成式(8)类型的高阶大惯性对象。
1.3.2 系统结构
分散式ADRC将回路之间的耦合视为扰动，并分别对每条回路单独设计ESO进行观测。具体地，忽略系统传递函数矩阵中的非对角传递函数，针对回路i对应的对角传递函数Gii进行ADRC控制器设计，对每条回路的控制器设计完成之后，即完成了对于多变量系统的分散式ADRC的设计。
对于系统传递函数矩阵中的对角传递函数为高阶大惯性环节的被控对象，考虑将补偿ADRC的设计方法应用在分散式ADRC当中。由于分散式ADRC各回路的设计是相互独立的，所以加入补偿环节后的控制器设计方法和参数整定过程还是与上文保持一致。称补偿ADRC与分散式ADRC结合后的控制方法为分散式补偿ADRC，以2×2阶多变量系统为例，其结构如图 3所示。图 3中，Gcp1、Gcp2分别为针对回路1和回路2设计的补偿环节；U =[U1 U2]T为系统的控制量；U ′=[U′ 1 U′ 2]T为补偿后的广义控制量。
1.3.3 系统传递函数推导
为对分散式补偿ADRC进行理论分析，现推导其系统闭环传递函数。以2×2阶多变量系统的1阶分散式补偿ADRC控制器为例，推导整个系统的闭环传递函数。先讨论补偿ADRC的结构，由式(5)可得ESO的状态方程：
\left[\begin{array}{c} \dot{z}_1 \\ \dot{z}_2 \end{array}\right]=\left[\begin{array}{cc} -\beta_1 & 1 \\ -\beta_2 & 0 \end{array}\right]\left[\begin{array}{l} z_1 \\ z_2 \end{array}\right]+\left[\begin{array}{cc} b_0 & \beta_1 \\ 0 & \beta_2 \end{array}\right]\left[\begin{array}{l} u_0 \\ y \end{array}\right] (12) 式(12)中，u0为图 2中的广义控制量。由于控制器为1阶，故式(9)中的m=1，将补偿环节式(9)加入状态反馈控制律式(6)中，可得：
u_0=\frac{k_{\mathrm{p}}\left(r-z_1\right)-z_2}{b_0} \frac{1}{(T s+1)^{n-1}} (13) 将式(13)代入式(12)中，可得：
\begin{aligned} {\left[\begin{array}{c} \dot{z}_1 \\ \dot{z}_2 \end{array}\right] } & =\left[\begin{array}{cc} -\beta_1-\frac{k_{\mathrm{p}}}{(T s+1)^{n-1}} & 1-\frac{1}{(T s+1)^{n-1}} \\ -\beta_2 & 0 \end{array}\right]\left[\begin{array}{l} z_1 \\ z_2 \end{array}\right]+ \\ & =\left[\begin{array}{cc} \frac{k_{\mathrm{p}}}{(T s+1)^{n-1}} & \beta_1 \\ 0 & \beta_2 \end{array}\right]\left[\begin{array}{l} r \\ y \end{array}\right] \\ & =\boldsymbol{A}_{\mathrm{c}}\left[\begin{array}{l} z_1 \\ z_2 \end{array}\right]+\boldsymbol{B}_{\mathrm{c}}\left[\begin{array}{l} r \\ y \end{array}\right] \end{aligned} (14) 其中，向量[z1 z2]T和[r y]T的系数矩阵依次定义为Ac和Bc。将[z1 z2]T变换至等式左边，得到向量[z1 z2]T与[r y]T之间的关系：
\left[\begin{array}{l} z_1 \\ z_2 \end{array}\right]=\left[s \boldsymbol{I}-\boldsymbol{A}_{\mathrm{c}}\right]^{-1} \boldsymbol{B}_{\mathrm{c}}\left[\begin{array}{l} r \\ y \end{array}\right] (15) 式(15)中，I为单位矩阵。利用线性系统的叠加原理，分别将式(15)中的r和y设置为0，则可得到z1和z2分别用y和r表示的两组形式。将其代入控制律u=(kp(r-z1)-z2)/b0，分别得到系统控制量u与系统输出y，系统控制量u和系统输入r之间的关系式：
\left\{\begin{array}{l} \frac{u(s)}{r(s)}=\frac{k_{\mathrm{p}}\left(s^2+\beta_1 s+\beta_2\right)}{b_0\left[s^2+\beta_1 s+k_{\mathrm{p}} s(T s+1)^{1-n}+\beta_2-\beta_2(T s+1)^{1-n}\right]} \\ \frac{u(s)}{y(s)}=-\frac{\left(k_{\mathrm{p}} \beta_1+\beta_2\right) s+k_p \beta_2}{b_0\left[s^2+\beta_1 s+k_{\mathrm{p}} s(T s+1)^{1-n}+\beta_2-\beta_2(T s+1)^{1-n}\right]} \end{array}\right. (16) 由图 4可知：
\left\{\begin{array}{l} G_{\mathrm{c}}(s)=\frac{u(s)}{r(s)} \\ F(s)=-\frac{r(s)}{y(s)}=-\frac{\frac{u(s)}{y(s)}}{\left(\frac{u(s)}{r(s)}\right)} \end{array}\right. (17) 将式(16)代入式(17)可得：
\left\{\begin{array}{l} G_c(s)=\frac{k_{\mathrm{p}}\left(s^2+\beta_1 s+\beta_2\right)}{b_0\left[s^2+\beta_1 s+k_{\mathrm{p}} s(T s+1)^{1-n}+\beta_2-\beta_2(T s+1)^{1-n}\right]} \\ F(s)=\frac{\left(k_{\mathrm{p}} \beta_1+\beta_2\right) s+k_{\mathrm{p}} \beta_2}{k_{\mathrm{p}}\left(s^2+\beta_1 s+\beta_2\right)} \end{array}\right. (18) 系统中的每条回路均可按上述运算过程转化为2自由度的形式，其中除了带宽参数不同外，各条回路中的传递函数Gc与F的形式均相同，而带宽参数是运用上文所述的补偿ADRC参数整定方法分别对各自回路的对角传递函数整定得到的。简化后的系统结构如图 5所示。
![]()
图 5 分散式补偿ADRC系统2自由度结构图Fig. 5. Two degree of freedom structure diagram of decentralized compensated ADRC system图 5中的输入向量R、控制向量U和输出向量Y均为2维向量，且有矩阵：
\boldsymbol{G}=\left[\begin{array}{ll} G_{11} & G_{12} \\ G_{21} & G_{22} \end{array}\right], \boldsymbol{G}_c=\left[\begin{array}{cc} G_{\mathrm{c} 1} & 0 \\ 0 & G_{\mathrm{c} 2} \end{array}\right], \boldsymbol{F}=\left[\begin{array}{cc} F_1 & 0 \\ 0 & F_2 \end{array}\right] (19) 其中，G为被控的2×2阶多变量系统矩阵，Gci和Fi表示将第i条回路中的控制器参数和对应的对角传递函数的参数代入式(18)计算得到的传递函数。根据图 5可以求得分散式补偿ADRC控制多变量系统的闭环传递函数Gd(s)为
\boldsymbol{G}_{\mathrm{d}}(s)=\frac{\boldsymbol{Y}(s)}{\boldsymbol{R}(s)}=\left(\boldsymbol{I}+\boldsymbol{G} \boldsymbol{G}_{\mathrm{c}} \boldsymbol{F}\right)^{-1} \boldsymbol{G} \boldsymbol{G}_{\mathrm{c}} (20) 2. 多变量系统的稳定区域判定
在多变量的频域分析方法中，逆奈奎斯特阵列设计方法得到了广泛的应用，而Gershgorin带和Ostrowski带是逆奈奎斯特阵列设计方法的基础。下面将分别介绍Gershgorin带和Ostrowski带，以及基于逆奈奎斯特阵列方法的稳定区域判定。
2.1 Gershgorin带与Ostrowski带
定义m×m阶复数矩阵Ap第i行的行估计值为
d_i=\sum\limits_{j=1}^m\left|a_{i j}\right|, j \neq i (21) 式中， \left|a_{i j}\right|表示矩阵Ap中第i行第j列的复数元的模。以矩阵Ap的对角元aii在复平面中的点为圆心， di为半径作圆，称为矩阵Ap第i行的行Gershgorin圆。
对于m×m阶传递函数矩阵Ap(s)，其内部每个元素都是随s变化的。当s随着频率由小到大逐渐变化时，Ap(s)的m个行Gershgorin圆也会随之发生变化，所以这些圆就会在复平面内扫出m条带状区域，称其为Ap(s)的行Gershgorin带。
若复数矩阵Ap为对角优势矩阵，定义Ap的第i行的行半径系数θi为
\left|a_{i j}\right|<\theta_i\left|a_{i i}\right| (22) 定义Ap的第i行的行压缩因子为φi，其值等于第i行外Ap的其他各行的行半径系数的最大值：
\boldsymbol{\varphi}_i=\max _{j(j \neq i)}\left\{\boldsymbol{\theta}_j\right\} (23) 与Gershgorin圆的不同之处在于，Ostrowski圆的半径是在Gershgorin圆半径的基础上乘以了该行的压缩因子，即矩阵的第i行的Ostrowski圆的半径 \hat{d}_i=\varphi_i d_i。故对角优势矩阵的Ostrowski带通常是在Gershgorin带的内部。
2.2 基于Ostrowski带的稳定区域判定
考虑图 6所示的多变量系统，其中，Q (s)为前向通道的传递函数矩阵，Fc=diag{f1，f2，…，fm}为反馈增益矩阵且各fi为非零实常数。
若开环系统稳定且Q (s)的逆矩阵 \hat{\boldsymbol{Q}}(s)每行的行Ostrowski带既不含每行对应的反馈增益矩阵Fc中的-fi点，也不包围-fi点，则闭环系统稳定[19]。
通常反馈增益矩阵Fc中fi数值的选择是在逆奈奎斯特阵列设计过程的最后一步：在设计完成整个前向通道后，画出 \hat{\boldsymbol{Q}}(s)的Ostrowski带，根据其与负实轴的交点中离原点较近的点到原点的距离确定fi的选择范围。根据上述稳定判据条件可知，该交点距离原点的位置越远，使系统稳定的反馈增益参数的选择区域就越大，即稳定区域越大。延伸这一思路，可以将此稳定区域大小的判定方法运用在对其他多变量系统设计方法的分析当中。
对于其他的多变量系统设计方法，系统设计的完成就意味着系统解耦的完成，即系统的闭环传递函数矩阵Gd(s)转化为了对角优势矩阵，故可以画出Gd(s)逆矩阵的Ostrowski带。依据不同系统逆矩阵的Ostrowski带就可以直观地对比使用不同方法设计系统的稳定区域大小，即可从稳定区域大小的角度定量分析不同设计方法的性能优劣。
3. 仿真实例
3.1 中速磨煤机多变量系统模型
考虑中速磨煤机系统模型[20]：
\left[\begin{array}{l} Y_1 \\ Y_2 \end{array}\right]=\left[\begin{array}{cc} \frac{1}{(20 s+1)^3} & \frac{1}{(25 s+1)^3} \\ \frac{1}{(80 s+1)^3} & \frac{-1}{(60 s+1)^3} \end{array}\right]\left[\begin{array}{l} R_1 \\ R_2 \end{array}\right] 分别使用分散式ADRC和分散式补偿ADRC对其进行控制，两种控制器的阶次均选择1阶。由于分散式ADRC本身就具有逆解耦能力，所以直接忽略被控对象的非对角传递函数，根据对角传递函数分别设定回路1(R1→Y1)和回路2(R2→Y2)的控制器参数。对于分散式补偿ADRC，设计回路1的补偿环节为Gcp1=1/(20s+1)2，设计回路2的补偿环节为Gcp2=1/(60s+1)2，再运用补偿ADRC的参数整定方法计算和调整得到控制器参数。两种方法的控制器参数如表 1所示。
表 1 控制系统参数Tab. 1. Control system parameters控制系统 回路 b0 ωc ωo 分散式ADRC 回路1 0.192 7 0.029 1 1.164 0 回路2 -0.064 2 0.009 7 0.038 8 分散式补偿ADRC 回路1 0.050 0 0.050 0 0.500 0 回路2 -0.031 3 0.016 7 0.167 0 依据上文推导的系统传递函数，分别绘制分散式ADRC与分散式补偿ADRC两种方法的系统逆矩阵Ostrowski带，如图 7所示。
![]()
图 7 分散式ADRC和分散式补偿ADRC的Ostrowski带Fig. 7. Ostrowski bands of decentralized ADRC and decentralized compensated ADRC以此Ostrowski带为据绘制出两个系统的稳定区域，如图 8所示。由图 8可以看出，分散式补偿ADRC的稳定区域大于分散式ADRC，说明分散式补偿ADRC对于高阶多变量系统有着更好的控制效果。
进行阶跃扰动实验，得到分散式ADRC和分散式补偿ADRC的控制效果如图 9所示。由图可以看出，相比于分散式ADRC，分散式补偿ADRC对于两条回路都具有更好的控制效果。
![]()
图 9 分散式ADRC和分散式补偿ADRC控制效果对比Fig. 9. Comparison of control effect between decentralized ADRC and decentralized compensated ADRC为分析在被控对象模型发生变化时两种方法的鲁棒性能，运用蒙特卡洛算法进行随机实验。在保持控制器参数不变的前提下，被控对象的每个传递函数中的时间常数T和增益K在90%~110%的范围内发生摄动，进行500次仿真实验，分别记录两条回路中设定值跟踪的超调量与调节时间。仿真结果如图 10所示。图中的点越靠近原点就拥有更好的动态性能，而同类型的点越集中表示该控制方法的鲁棒性能更强。由图中结果可知，在被控对象的参数发生摄动时，分散式补偿ADRC的性能散点图分布范围更小，并更靠近坐标轴右下方，即分散式补偿ADRC拥有更好的动态性能和鲁棒性。
![]()
图 10 分散式ADRC与分散式补偿ADRC蒙特卡洛实验图Fig. 10. Monte Carlo experiment diagram of decentralized ADRC and decentralized compensated ADRC3.2 气体流量装置系统模型
考虑气体流量装置实验管路流量、压力耦合系统模型[21]。系统的输入为两个调节阀的给定电流I1和I2，输出为管道的流量Qm和压力p。该系统为一个2输入2输出系统：
\begin{aligned} & {\left[\begin{array}{l} Q_{\mathrm{m}} \\ p \end{array}\right] } \\ = & {\left[\begin{array}{ll} G_{Q_m I_1} & G_{Q_m I_2} \\ G_{p I_1} & G_{p l_2} \end{array}\right]\left[\begin{array}{l} I_1 \\ I_2 \end{array}\right] } \\ = & {\left[\begin{array}{cc} \frac{48.252}{0.36 s^2+0.9163 s+1} & \frac{55.913}{0.348 s^2+0.905 s+1} \\ \frac{54.55}{0.0894 s^2+0.4258 s+1} & \frac{-70.429}{0.084 s^2+0.5254 s+1} \end{array}\right]\left[\begin{array}{l} I_1 \\ I_2 \end{array}\right] } \end{aligned} 为对其进行控制器设计，运用文[17]的近似方法，将系统传递函数矩阵中的主对角传递函数近似为高阶大惯性环节：
G_{Q_m I_1} \approx \frac{48.252}{(0.336 s+1)^3}, G_{\mathrm{P} l_2} \approx \frac{-70.429}{(0.178 s+1)^3} 故1阶分散式补偿ADRC的两个补偿环节分别为Gcp1=1/(0. 336s+1)2和Gcp2=1/(0. 178s+1)2。根据参数整定方法得到两个控制器参数如表 2所示。
表 2 不同控制器的控制参数Tab. 2. Control parameters of different controllers控制系统 回路 b0 ωc ωo 分散式ADRC 回路1 553 1.73 6.92 回路2 1 525 3.27 13.1 分散式补偿ADRC 回路1 209 2.98 29.8 回路2 -272 5.62 56.2 系统设计完成后，分别绘制两种方法的系统传递函数逆矩阵Ostrowski，如图 11所示。
由图 11可知，相较于分散式ADRC，分散式补偿ADRC对角传递函数的奈奎斯特曲线都有所抬高，且Ostrowski圆都有所缩小，这就使得其Ostrowski带与负实轴的交点往左移动，说明分散式补偿ADRC的对角优势更加明显且拥有更大的稳定区域。分别绘制两个系统的稳定区域如图 12所示，其中分散式补偿ADRC在两条回路中都拥有更大的稳定区域。
进行阶跃扰动实验，得到两种控制方法的控制效果如图 13所示。可以看到，分散式补偿ADRC的响应更快，并且其对于回路之间耦合的抑制作用也要强于分散式ADRC。
进行蒙特卡洛随机实验，保持各控制器的参数不变，使气体流量模型矩阵中的各个传递函数中的参数在90%~110%的范围内发生摄动，得到仿真结果如图 14所示。可以发现分散式补偿ADRC在被控模型发生扰动的时候拥有更好的动态性能。且在两条回路当中，分散式补偿ADRC的各点相较于分散式ADRC的各点均更加集中，说明分散式补偿ADRC有着更强的鲁棒性。
4. 结论
本文针对所提出的分散式补偿ADRC，通过理论分析得到了该方法控制多变量系统的闭环传递函数，并将多变量控制系统设计方法中的逆奈奎斯特阵列法应用在了系统的稳定区域大小的定量分析当中。通过实例仿真发现，在加入补偿环节后，分散式补偿ADRC在使用低阶控制器的情况下对高阶多变量系统有着更好的控制效果，其动态性能和鲁棒性能均优于分散式ADRC。
-
[1] Margarete A S,Ken P. Meeting the closed-loop challenge: The case of remanufacturing[J]. California Management Review,2004,46(2): 74-89. [2] Ferrer G,Swaminathan J M. Managing new and differentiated remanufactured products[J]. European Journal of Operational Research,2010,203(2): 370-379.
[3] Chen J M,Chang C I. The competitive strategy of a closed-loop supply chain with remanufacturing[J]. Transportation Research: Part E,2012,48(2): 387-400.
[4] 郭军华,李帮义,倪明. WTP差异下的再制造进入决策研究[J]. 中国管理科学,2013,21(1): 149-156. Guo J H,Li B Y,Ni M. Research on remanufacturing entry decision based on WTP differentiation[J]. Chinese Journal of Management Science,2013,21(1): 149-156.
[5] Giuntini R,Gaudette K. Remanufacturing: The next great opportunity for boosting US productivity[J]. Business Horizons,2003,46(6): 41-48.
[6] Atasu A,Jr Guide V D R,Van Wassenhove L N. So what if remanufacturing cannibalizes my new product sales?[J]. California Management Review,2010,52(2): 56-76.
[7] Martin P,Jr Guide V D R,Craighead C W. Supply chain sourcing remanufacturing operations: An empirical investigation of remake versus buy[J]. Decision Sciences,2010,41(2): 301-324.
[8] Savaskan R C,Bhattacharya S,Van Wassenhove L N. Closed-loop supply chain models with product remanufacturing[J]. Management Science,2004,50(2): 239-252.
[9] Debo L G,Toktay L B,Van Wassenhove L N. Market segmentation and product technology selection for remanufacturable products[J]. Management Science,2005,51(8): 1193-1205.
[10] 孙浩,达庆利. 基于产品差异的再制造闭环供应链定价与协调研究[J]. 管理学报,2010,7(5): 733-738. Sun H,Da Q L. Pricing and coordination of remanufacturing closed-loop supply chain based on product differentiation[J]. Chinese Journal of Management,2010,7(5): 733-738.
[11] 熊中楷,王凯,熊榆,等. 考虑经销商从事再制造的闭环供应链模式[J]. 系统工程学报,2011,26(6): 792-799. Xiong Z K,Wang K,Xiong Y,et al. Closed-loop supply chain considering that the distributor engages in remanufacturing[J]. Journal of Systems Engineering,2011,26(6): 792-799.
[12] Geyer R. The economics of remanufacturing under limited component durability and finite product life cycles[J]. Management Science,2007,53(1): 88-100.
[13] 李新军,达庆利. 基于产品寿命次数的闭环供应链生态优化问题研究[J]. 管理工程学报,2009,23(2): 50-54. Li X J,Da Q L. Monopoly production optimization problem for product lifetime in closed-loop supply chains[J]. Journal of Industrial Engineering and Engineering Management,2009,23(2): 50-54.
[14] Nagurney A,Dong J,Zhang D. A supply chain network equilibrium model[J]. Transportation Research: Part E,2002,38(5): 281-303.
[15] Dong J,Zhang D,Nagurney A. A supply chain network equilibrium model with random demands[J]. European Journal of Operational Research,2004,156(1): 194-212.
[16] 胡劲松,徐元吉,刘芳霞,等. 具有模糊需求的多商品流供应链网络均衡研究[J]. 控制与决策,2012,27(5): 665-672. Hu J S,Xu Y J,Liu F X,et al. A multi-products flow supply chain network equilibrium with fuzzy demand[J]. Control and Decision,2012,27(5): 665-672.
[17] Nagurney A,Cruz J,Dong J,et al. Supply chain networks,electronic commerce,and supply side and demand side risk[J]. European Journal of Operational Research,2005,164(1): 120-142.
[18] Hammond D,Beullens P. Closed-loop supply chain network equilibrium under legislation[J]. European Journal of Operational Research,2007,183(2): 895-908.
[19] Yang G F,Wang Z P,Li X Q. The optimization of the closed-loop supply chain network[J]. Transportation Research: Part E,2009,45(1): 16-28.
[20] 杨玉香,周根贵. EPR下供应链网络报废产品排放内生污染税模型[J]. 管理科学学报,2011,14(10): 67-76. Yang Y X,Zhou G G. An endogenous EOL products emission taxes model for supply chain network under EPR[J]. Journal of Management Sciences in China,2011,14(10): 67-76.
[21] 张桂涛,胡劲松,孙浩,等. 具有缺陷产品的双渠道闭环供应链网络均衡[J]. 中国管理科学,2013,21(5): 68-79. Zhang G T,Hu J S,Sun H,et al. The equilibrium of closed-loop supply chain network with dual channel and imperfect product [J]. Chinese Journal of Management Science,2013,21(5): 68-79.
[22] 张铁柱,周倩. 双渠道多期供应链网络均衡模型研究[J]. 计算机集成制造系统,2008,14(8): 1512-1520. Zhang T Z,Zhou Q. Multi-period supply chain network equilibrium model with two-channel[J]. Computer Integrated Manufacturing Systems,2008,14(8): 1512-1520.
[23] Hamdouch Y. Multi-period supply chain network equilibrium with capacity constraints and purchasing strategies[J]. Transportation Research: Part C,2011,19(5): 803-820.
[24] Korpelevich G M. The extra gradient method for finding saddle points and other problems[J]. Matekon,1977,13(1): 35-49.
[25] He B S,Xu Y,Yuan X M. A logarithmic-quadratic proximal prediction-correction method for structured monotone variational inequalities[J]. Computational Optimization and Applications,2006,35(1): 19-46.
下载:
计量
- 文章访问数: 874
- HTML全文浏览量: 4
- PDF下载量: 227
