Type-2 Adaptive Fuzzy Modeling and Oxygen Excess Ratio Control for PEMFC Air Supply System
-
摘要: 质子交换膜燃料电池（Proton exchange membrane fuel cell，PEMFC）空气供给系统存在外部扰动和参数不确定等动态特性，难以实现精准建模和控制.本文结合精确线性化和二型模糊逻辑系统，提出一种自适应控制器实现PEMFC空气供给系统的建模与过氧比控制.该控制器不需要PEMFC空气供给系统模型结构和参数完全已知的条件，而是通过二型模糊逻辑系统在线逼近PEMFC空气供给系统中的未建模动态并从Lyapunov函数中导出自适应参数，从而保证系统收敛性与稳定性.通过稳定性分析证明了该控制器作用下系统跟踪误差的有界性，仿真实验进一步验证了该控制器的有效性与实用性.
-
关键词:
- 二型模糊逻辑系统 /
- 自适应控制 /
- 精确线性化 /
- Lyapunov稳定性 /
- 过氧比
Abstract: Proton exchange membrane fuel cell (PEMFC) air supply system has the characteristics of external disturbances and uncertain parameters, which make it difficult to achieve accurate modeling and stability control. In this paper, an adaptive controller is proposed to control the oxygen excess ratio of PEMFC air supply system by using the type-2 fuzzy logic systems. The controller does not need the known conditions of PEMFC system model but approximates unmodeled dynamics in the system by the adaptive fuzzy system whose the parameter adjustment is derived based on the Lyapunov theory. The stability analysis shows that the system is stable under the control of the controller. Simulation results demonstrate the usefulness and effectiveness of our proposed control strategy. -
为了解决能源危机和全球变暖问题, 当前世界各国都在加大新能源开发力度.质子交换膜燃料电池(Proton exchange membrane fuel cell, PEMFC)是一种通过电化学反应将储存在燃料中的化学能直接转换为电能的电化学装置, 具有高功率、高能量转换效率、低温启动、无环境污染等优点, 应用前景广阔、已逐渐在新能源汽车和电力系统等领域得到实际应用[1-4].
合理控制PEMFC空气供给系统阴极注入空气流量与电化学反应消耗空气流量之比, 即过氧比(Oxygen excess ratio, OER)是提高PEMFC系统输出净功率的有效方法.为此, 国内外学者对PEMFC空气供给系统的建模与控制问题作了大量研究.文献[5-7]基于以下理想假设建立了PEMFC空气供给系统的基础数学模型: 1)视空气为理想气体, 遵循理想气体定律, 且干燥空气中$21\, \%$为氧气成分, $79\, \%$为氮气成分; 2)空压机机械效率较高, 摩擦力矩忽略不计; 3)冷却器体积忽略不计, 并可迅速降低空气温度到80℃; 4)排出阴极的空气特性与阴极内空气特性相同; 5)阴极内水蒸气平衡得到很好地控制, 且无液态水存在.基于以上基础模型及改进模型, 国内外学者提出采用线性反馈控制、模型预测控制和滑模控制等方法对PEMFC空气供给系统进行了控制研究[8-12].但以上文献控制器的设计均是基于理想模型, 且在模型结构、参数或者模型函数界严格已知条件下进行的, 目前如何有效辨识PEMFC模型结构、参数和函数界仍是一困难问题.实践表明, 模糊逻辑系统和自适应控制技术相结合是解决复杂而无法建立精确数学模型系统的控制问题的有效途径之一.但模糊逻辑系统和自适应控制技术在处理PEMFC空气供给控制系统中不确定性和不精确性问题的研究很少, 特别是采用二型模糊逻辑系统、自适应技术并结合Lyapunov稳定性分析开展研究的鲜有报道.
自从Zadeh教授提出Ⅰ型模糊集合(T1-fuzzy set, T1-FS)理论以来, 基于T1-FS的Ⅰ型模糊逻辑系统(T1-fuzzy logic system, T1-FLS)成为处理不确定性、模糊性和不精确性建模问题的有效方法, 并在众多领域得到了广泛应用.随着理论和实际应用的发展, 人们逐步意识到一型模糊集合在描述多重模糊不确定性方面的局限性, 使得一些从事一型模糊集合研究的学者们开始转而研究二型模糊集合.特别是近年来二型模糊集合和二型模糊系统的理论与应用得到了快速发展.这是因为一型模糊系统采用了由精确隶属度函数表示模糊集合, 其直接处理模糊规则不确定性的能力非常有限, 而二型模糊集合和二型模糊系统能较好地解决这些不确定性问题[13-17].虽然该领域在Mendel教授的推动下得到了快速发展, 但在许多原创性问题上还在被探索研究, 如最近两年我国学者王立新教授探索了二型模糊集合和二型模糊逻辑系统的新框架[18]、王飞跃和莫红教授探索了二型模糊集合的一系列新定义[19].这些新框架新定义与Mendel定义的流行框架哪个更好, 或者说针对不同应用对象哪个定义更匹配、目前没有准确对比报道还需进一步研究.
模糊系统和自适应控制的结合是对非线性不确定系统进行建模与控制的一种有效方法, 自适应模糊建模与控制的基本出发点是仿人的智能以实现对复杂不确定性系统进行有效的建模与控制, 它具有从环境自学习、适应环境的能力.该方法不需要对象模型结构和参数完全已知的条件, 仅根据被控系统结构关系即可实现有效控制, 已在许多领域的建模和控制中得到广泛关注[20-25].目前, 解决PEMFC空气供给系统的控制问题, 一般都是首先建立空气供给系统的机理模型、然后进行控制器设计.而进行智能建模和控制方面的研究工作较少.
本文对PEMFC空气供给系统的模糊建模与过氧比自适应控制问题做了一些探索研究、主要工作如下: 1)本文首先对国际文献中的一些推导错误及多漏写失误进行了改正, 从而保证仿真状态空间模型的正确性; 2)为了设计非模型的自适应模糊控制器, 采用精确线性化技术求解了系统相对阶数; 3)基于等价模型, 对PEMFC空气供给系统模型中的多个未建模动态采用二型模糊系统建模, 这种方法对PEMFC空气供给系统进行建模在国内外很少有报道; 4)结合具体的PEMFC空气供给系统的建模与控制问题, 从Lyapunov稳定性导出自适应律并设计了相应自适应控制器.这种结合二型模糊系统建模、自适应控制和稳定性分析, 在国际主要文献上未见.大部分PEMFC系统文献局限在简单的模糊控制、或者有自适应功能但缺少稳定性分析等; 5)仿真及其对比实验验证了该控制器的有效性与实用性.
1. PEMFC空气供给系统的数学模型
PEMFC系统主要由氢气供给系统、空气供给系统和反应电堆三部分构成(见图 1), 其中PEMFC空气供给系统包括空压机、供应管道、冷却器、加湿器、电堆阴极以及回流管道等组成. PEMFC空气供给系统具体结构如图 1虚框所示.
1.1 PEMFC空气供给系统的建模
为了实现PEMFC空气供给系统的高级控制方法的研究.文献[5]在理想假设基础上, 基于文献[7]的流体力学、热力学和电化学等基础物理方程建立了含6个状态变量的PEMFC空气供给系统的非线性状态空间模型, 即:
$ \begin{equation} \begin{cases} \dot{x}_{1} =\frac{1}{J_{cp}}[\tau_{cm}(u, x_{1})-\tau_{cp}(x_{1}, x_{2})]\\ \dot{x}_{2} = \frac{\gamma R_{a}}{V_{sm}}[T_{cp}(x_{2})W_{cp}(x_{1}, x_{2})-\\ \ \ \ \ \ \ \ T_{sm}(x_{2}, x_{3})W_{sm, \rm out}(x_{2}, x_{4}, x_{5})]\\ \dot{x}_{3} =W_{cp}(x_{1}, x_{2})-W_{sm, \rm out}(x_{2}, x_{4}, x_{5})\\ \dot{x}_{4} =W_{{\rm O_{2}}, \rm in}(x_{2}, x_{4}, x_{5})-W_{{\rm O_{2}}, \rm react}({\rm I}_{st})-\\ \ \ \ \ \ \ \ W_{{\rm O_{2}}, \rm out}(x_{4}, x_{5}, x_{6})\\ \dot{x}_{5} =W_{{\rm N_{2}}, \rm in}(x_{2}, x_{4}, x_{5})-W_{{\rm N_{2}}, \rm out}(x_{4}, x_{5}, x_{6}) \\ \dot{x}_{6} =\frac{R_{a}T_{fc}}{V_{rm}M_{a}}[ W_{ca, \rm out}(x_{4}, x_{5}, x_{6})-W_{rm, \rm out}(x_{6})] \end{cases} \end{equation} $
(1) 其中, $J_{cp}$为空压机转动惯量, $u$为空压机工作电压, $\tau_{cm} $为空压机驱动力矩, $\tau_{cp}$为空压机负载力矩, $\gamma$为空气比热容比, $R_{a}$为空气气体常数, $V_{sm}$为供给管道体积, $T_{cp}$为空压机流出空气温度, $W_{cp}$为空压机流出空气流量, $T_{sm}$为供给管道流出空气温度, $W_{sm, out}$为供给管道流出空气流量, $W_{{\rm O_{2}}, \rm{in}}$为阴极注入氧气流量, $W_{{\rm O_{2}}, \rm{react}}({\rm I}_{st})$为电堆反应消耗氧气流量, ${\rm I}_{st}$为负载电流, $W_{{\rm O_{2}}, \rm{out}}$为阴极流出氧气流量, $W_{{\rm N_{2}}, \rm{in}}$为阴极注入氮气流量, $W_{{\rm N_{2}}, \rm{out}}$为阴极流出氮气流量, $T_{fc}$为电堆温度, $V_{rm}$为回流管道体积, $W_{ca, \rm{out}}$为阴极流出空气流量, $W_{rm, \rm{out}}$为回流管道流出空气流量, $M_{a}$为空气摩尔质量.
注1. 该模型中存在多处推导失误公式、参数及误打的错误需要修正, 详细修正内容见本文附录A.
基于方程组(1)和附录A的修正后, PEMFC空气供给系统可表示为如下状态空间方程:
$ \begin{align} \dot{\boldsymbol x}={\boldsymbol f}(\boldsymbol x)+{\boldsymbol g}\cdot u+{\boldsymbol d}\cdot {\rm I}_{st} \end{align} $
(2) 其中, $\boldsymbol x$表示该系统的状态向量、其含义如表 1所示, PEMFC系统负载电流${\rm I}_{st}$视为外界干扰, ${\boldsymbol f}(\boldsymbol x)$、${\boldsymbol g}$和${\boldsymbol d}$分别为:
表 1 PEMFC空气供给系统状态变量Table 1 State variables of PEMFC air supply system状态变量 符号 单位 空压机转速 $x_{1}=\omega_{cp}$ rad/s 供给管道内空气压强 $x_{2}=P_{sm}$ Pa 供给管道内空气质量 $x_{3}=m_{sm}$ kg 阴极内氧气质量 $x_{4}=m_{{\rm O_{2}}}$ kg 阴极内氮气质量 $x_{5}=m_{{\rm N_{2}}}$ kg 回流管道内空气压强 $x_{6}=P_{rm}$ Pa $ {\boldsymbol f}(\boldsymbol x)=\left[ \begin{array}{cc} f_{1}(x_{1}, x_{2}) \\ f_{2}(x_{1}, x_{2}, x_{3}, x_{4}, x_{5})\\ f_{3}(x_{1}, x_{2}, x_{4}, x_{5})\\ f_{4}(x_{2}, x_{4}, x_{5}, x_{6})\\ f_{5}(x_{2}, x_{4}, x_{5}, x_{6})\\ f_{6}(x_{4}, x_{5}, x_{6})\\ \end{array} \right] $
$ {\boldsymbol g}=\left[\begin{array}{cccc} \frac{\eta_{cm}K_{t}}{J_{cp}R_{cm}} \\ 0 \\ 0 \\ 0 \\ 0 \\ 0 \\ \end{array}\right], {\boldsymbol d}=\left[\begin{array}{cccc} 0 \\ 0 \\ 0 \\ -\frac{nM_{{\rm O_{2}}}}{4F} \\ 0 \\ 0 \\ \end{array}\right] $
其中, $\eta_{cm}$为空压机电机的机械效率, $K_{t}$、$R_{cm}$为空压机电机参数, $n$为单体电池个数, $M_{{\rm O_{2}}}$为氧气摩尔质量, $F$为法拉第常数. $f_{1}(x_{1}, x_{2})-f_{6}(x_{4}, x_{5}, x_{6})$函数的具体表达式见附录B.
根据近年来国内外研究可知, PEMFC系统负载电流一定时, 当系统过氧比$\lambda_{{\rm O_{2}}}$达到理想过氧比$\lambda_{{\rm O_{2}}, ref}$时, 系统输出净功率为最大值(见图 2).因此, 选择系统输出方程为:
图 2 PEMFC系统负载电流、过氧比和输出净功率的关系Fig. 2 Power relationship of PEMFC system load current, OER, and output net$ \begin{equation} \begin{cases} &y=h(\boldsymbol x)=\lambda_{{\rm O_{2}}}\\ &\lambda_{{\rm O_{2}}}=\frac{W_{{\rm O_{2}}, \rm in}(x_{2}, x_{4}, x_{5})}{W_{{\rm O_{2}}, \rm react}({\rm I}_{st})} \end{cases} \end{equation} $
(3) PEMFC系统理想过氧比$\lambda_{{\rm O_{2}, ref}}$可表示为[8]:
$ \begin{align} y_{m}=\lambda_{{\rm O_{2}, ref }}& = 5\times10^{-8 }{\rm I}_{st}^{3}-2.87\times10^{-5 }{\rm I}_{st}^{2}+ \notag\\ & 2.23\times10^{-3 }{\rm I}_{st}+2.5 \end{align} $
(4) 结合式(2)和式(3), PEMFC空气供给系统状态空间表达式可表示为:
$ \begin{align} \left\{\begin{aligned} &\dot{\boldsymbol x}={\boldsymbol f}(\boldsymbol x)+{\boldsymbol g}\cdot u+{\boldsymbol d}\cdot {\rm I}_{st}\\ & y=h(\boldsymbol x) \end{aligned}\right. \end{align} $
(5) 1.2 PEMFC空气供给系统的精确线性化
本文通过非线性状态和反馈变换, 对式(5)的非线性状态实现精确线性化, 输出方程李导数$L_{{\boldsymbol f}}h(\boldsymbol{x})$、$L_{{\boldsymbol g}}h(\boldsymbol{x})$和$L_{{\boldsymbol d}}h(\boldsymbol{x})$可表示为:
$ \begin{align} \left\{\begin{aligned} & L_{{\boldsymbol f}}h(\boldsymbol{x})=f_{2}\frac{\partial h(\boldsymbol x)}{\partial x_{2}}+f_{4}\frac{\partial h(\boldsymbol x)}{\partial x_{4}}+f_{5}\frac{\partial h(\boldsymbol x)}{\partial x_{5}}\\ & L_{{\boldsymbol g}}h(\boldsymbol{x})=0\\ & L_{{\boldsymbol d}}h(\boldsymbol{x})=-\frac{nM_{\rm O_{2}}}{4F} \frac{\partial h(\boldsymbol x)}{\partial x_{4}} \end{aligned}\right. \end{align} $
(6) 由于$L_{{\boldsymbol g}}h(\boldsymbol{x})=0$, 对系统继续求解李导数. $L^{2}_{{\boldsymbol f}}h(\boldsymbol{x})$、$L_{{\boldsymbol g}}L_{{\boldsymbol f}}h(\boldsymbol{x})$和$L_{{\boldsymbol d}}L_{{\boldsymbol f}}h(\boldsymbol{x})$分别表示为:
$ \begin{align} \left\{\begin{aligned} & L^{2}_{{\boldsymbol f}}h(\boldsymbol{x}) =f_{1}\frac{\partial L_{{\boldsymbol f}}h(\boldsymbol{x})}{\partial x_{1}}+f_{2}\frac{\partial L_{{\boldsymbol f}}h(\boldsymbol{x})}{\partial x_{2}} + \\ & \ \ \ \ \ \ \ \ \ f_{3}\frac{\partial L_{{\boldsymbol f}}h(\boldsymbol{x})}{\partial x_{3}} + f_{4}\frac{\partial L_{{\boldsymbol f}}h(\boldsymbol{x})}{\partial x_{4}}+\\ & \ \ \ \ \ \ \ \ \ f_{5}\frac{\partial L_{{\boldsymbol f}}h(\boldsymbol{x})}{\partial x_{5}}+f_{6}\frac{\partial L_{{\boldsymbol f}}h(x)}{\partial x_{6}}\\ & L_{{\boldsymbol g}}L_{{\boldsymbol f}}h(\boldsymbol{x}) =\frac{\eta_{cm} K_{t}}{J_{cp} R_{cm}} \frac{\partial L_{{\boldsymbol f}}h(\boldsymbol{x})}{\partial x_{1}}\neq 0\\ & L_{{\boldsymbol d}}L_{{\boldsymbol f}}h(\boldsymbol{x}) = -\frac{nM_{\rm O_{2}}}{4F} \frac{\partial L_{{\boldsymbol f}}h(\boldsymbol{x})}{\partial x_{4}} \end{aligned}\right. \end{align} $
(7) 由上可知PEMFC空气供给系统相对阶数为2, 小于表达式(5)所描述的系统阶数, 方程(5)可转化成以下等价系统:
$ \begin{align} \left\{\begin{aligned} & \dot{z}_{1}=z_{2}\\ & \dot{z}_{2}=\rho(\boldsymbol{x})+\varphi(\boldsymbol{x})u+\vartheta(\boldsymbol{\chi}) \end{aligned}\right. \end{align} $
(8) 其中, $\boldsymbol{\boldsymbol{\chi}}=[\boldsymbol{x}, {\rm I}_{st}]^{{\rm{T}}}$, $ \rho(\boldsymbol{x})=L^{2}_{{\boldsymbol f}}h(\boldsymbol{x}) $, $\varphi(\boldsymbol{x})= L_{{\boldsymbol g}}L_{{\boldsymbol f}}h(\boldsymbol{x})$, $\vartheta(\boldsymbol{\chi})=L_{{\boldsymbol d}}L_{{\boldsymbol f}}h(\boldsymbol{x}){\rm I}_{st}+\frac{\rm{d}}{{\rm d}t}(L_{{\boldsymbol d}}h(\boldsymbol{x}){\rm I}_{st})$.
对于系统(8), 如果$\rho(\boldsymbol{x})$、$\varphi(\boldsymbol{x})$和$\vartheta(\boldsymbol{\chi})$结构和参数完全已知, 则理想控制律为:
$ \begin{align} \left\{\begin{aligned} & u=\varphi(\boldsymbol{x})^{-1}(v-\rho(\boldsymbol{x})-\vartheta(\boldsymbol{\chi}))\\ & v=\ddot{y}_{m}+k_{1}\dot{{ e}}+k_{2}{e} \end{aligned}\right. \end{align} $
(9) 其中, $e=y_{m}-y$.然而在$\rho(\boldsymbol{x})$、$\varphi(\boldsymbol{x})$和$\vartheta(\boldsymbol{\chi})$未知情况下, 获得理想控制器(9)是不可能的.在这种情况下, 本文设计下面Ⅱ型模糊逻辑系统分别逼近未建模动态$\rho(\boldsymbol{x})$、$\varphi(\boldsymbol{x})$和$\vartheta(\boldsymbol{\chi})$.
2. Ⅱ型模糊逻辑系统
只要在模糊规则前件或后件隶属度函数中包含二型模糊集(T2-FS), 相应的模糊系统即称为二型模糊逻辑系统(T2-FLS).本文模糊规则前件采用二型模糊集(T2-FS), 而后件由于采用自适应调节、选用一型模糊集(T1-FS), 该二型模糊系统的降型采用Karnik和Mendel提出的方法(国际上简称KM算法[26]).由于本文所采用的二型模糊系统后件是一型自适应模糊集(T1-FS), 故在二型模糊系统的降型处理时仅调用一次KM算法、从而提高了算法在线实时计算效率.
定义1[27-29].一个二型模糊集合(见图 3) $\tilde{A}$可以表示为一个二型隶属函数$\tilde{A}(x, u)$, 其中, $x\in X $, $u\in J_{x}\subseteq[0, 1]$, 即:
$ \begin{equation} \tilde{A}=\{((x, u), \mu_{\tilde{A}}(x, u))|x\in {X}, u\in J_{x}\subseteq[0, 1]\} \end{equation} $
(10) 也可表示下面形式:
$ \begin{equation} \tilde{A}=\int_{x\in {X}}\int_{u\in J_{x}}\frac{\mu_{\tilde{A}}(x, u)}{(x, u)}{\rm d}u{\rm d}x \end{equation} $
(11) 其中, $x$为主变量, $u$为次变量, $J_{x}\subseteq[0, 1]$为主隶属度, $0\leq\mu_{\tilde{A}}(x, u)\leq1$为次隶属度, $\int\int$表示所有可允许$x$与$u$之并.由于本文采用的是区间二型模糊集, 取$\mu_{\tilde{A}}(x, u)=1.$针对上面二型模糊集合$\tilde{A}$的定义, 其不确定覆盖域(Footprint of uncertainty, FOU)定义如下:
$ \begin{equation} {\rm{FOU}}(\tilde{A})=\bigcup\limits_{x\in {X}} J_{x} \end{equation} $
(12) 注2. 上面Mendel对二型模糊集合和不确定覆盖域的定义, 被国外绝大部分文献采用.但我国学者莫红和王飞跃教授认为此定义有欠缺, 对二型模糊集合$\tilde{A}$和不确定覆盖域(FOU)进行了新定义.文献[30]认为一个二型模糊集合$\tilde{A}$的不确定覆盖域(FOU)为论域上每一点与在该点的主隶属度的笛卡尔积之并, 也就是:
$ \begin{equation} {\rm{FOU}}(\tilde{A})=\bigcup\limits_{x\in {X}} x\times J_{x}. \end{equation} $
(13) 如果从数学角度看, 式(13)比式(12)更严谨.假设按照式(12)计算FOU, 结果可能是一个区间, 而非Mendel所想表达的图 3阴影面积.造成这种的原因可能是Mendel先预设定了不确定覆盖域(FOU), 然后用符号$\bigcup$表示累加移动的意思, 而非传统意义上的并运算.本文不确定覆盖域(FOU)的计算按照式(13), 或者用式(12)计算、但符号$\bigcup$表示累加移动的意思.
在上面二型模糊集合定义后, 考虑$n$输入单输出的Mamdani区间二型模糊逻辑系统, 其中输入$x_{i}\in X_{i}$, 输出$y\in Y$, $i=1, \cdots, n$.假设IT2-FLS中包含如下$M$条模糊规则:
$ \begin{eqnarray} R^{(j)}: {\rm{IF}}\ x_{1}\ \ {\rm is} \ \tilde{A}_{1}^{j}, \ \cdots, x_{n}\ \ {\rm is}\ \ \tilde{A}_{n}^{j}, \ \nonumber\\ \ {\rm{THEN}}\ \ y \ {\rm is} \ B^{j}, j=1, \cdots, M \end{eqnarray} $
(14) 其中, $\tilde{A}_{i}^{j}$与~$B^{j}$分别表示前提变量~$x_{i}$与结论变量$y$在第$j$条规则所对应的二型和一型模糊集合.对于上面的Ⅱ型模糊系统(14), 基于KM算法的降型在线计算过程如下:
步骤1. 对于某条规则, 计算输入向量$\boldsymbol{x}$的激活区间:
$ \begin{equation} \tilde{A}^{j}:\left\{ \begin{array}{lll} A^{j}(\boldsymbol{x})\equiv[\underline{a}^{j}(\boldsymbol{x}), \overline{a}^{j}(\boldsymbol{x})]\\ \underline{a}^{j}(\boldsymbol{x})\equiv T_{i=1}^{n}\underline{\mu}_{\tilde{A}_{i}^{j}}(x_{i}) \\ \overline{a}^{j}(\boldsymbol{x})\equiv T_{i=1}^{n}\overline{\mu}_{\tilde{A}_{i}^{j}}(x_{i}) \end{array} \right. \end{equation} $
(15) 步骤2. 对于某条规则, 输出变量$y$的最大隶属度所对应的值为:
$ \begin{equation} C_{B^{j}}=\bar{y}_{j}=\Theta_{j} \end{equation} $
(16) 步骤3. 在前件的激活区间$\tilde{A}^{j}$和后件的输出$\Theta_{j}$计算后, 基于KM算法的降型可表示如下:
$ \begin{equation} Y_{ACOS}:\left\{ \begin{array}{lll} Y_{ACOS}(\boldsymbol{x})=\frac{1}{[y_{l}(\boldsymbol{x}), y_{r}(\boldsymbol{x})]}\\ y_{l}(\boldsymbol{x})=\frac{\min\limits_{ \omega^{l}_{j}\in [\underline{a}^{j}(\boldsymbol{x}), \overline{a}^{j}(\boldsymbol{x})]}\sum\limits_{j=1}^{ M}\Theta_{j}\omega^{l}_{j}}{\sum\limits_{j=1}^{ M}\omega^{l}_{j}} \\ y_{r}(\boldsymbol{x})=\frac{\max\limits_{ \omega^{r}_{j}\in [\underline{a}^{j}(\boldsymbol{x}), \overline{a}^{j}(\boldsymbol{x})]}\sum\limits_{j=1}^{ M}\Theta_{j}\omega^{r}_{j}}{\sum\limits_{j=1}^{ M}\omega^{r}_{j}} \end{array} \right. \end{equation} $
(17) 其中, $[y_{l}, y_{r}]$为用KM算法得到的输出区间.
步骤4. 加权平均反模糊化后, 可将上面的IT2-FLS简化成下面基函数的形式:
$ \begin{equation} y=\frac{1}{2}\boldsymbol{\Theta}^{{\rm{T}}}[\boldsymbol{\xi}_{l}(\boldsymbol{x})+\boldsymbol{\xi}_{r}(\boldsymbol{x})] \end{equation} $
(18) 其中, ${\boldsymbol\Theta}=[{\Theta}_{1}, \cdots, {\Theta}_{M}]^{{\rm{T}}}$为二型模糊系统的自适应参数向量, 二型模糊系统的左基函数为$\boldsymbol{\xi}_{l}(\boldsymbol{x})= [\omega^{l}_{1}/\sum\nolimits_{j=1}^{ M}\omega^{l}_{j}, \cdots, \omega^{l}_{M}/\sum\nolimits_{j=1}^M\omega^{l}_{j}]^{{\rm{T}}}$, 二型模糊系统的右基函数为$\boldsymbol{\xi}_{r}(\boldsymbol{x})=[\omega^{r}_{1}/\sum\nolimits_{j=1}^{ M}\omega^{r}_{j}, \cdots, \omega^{r}_{M}/\sum\nolimits_{j=1}^{ M}\omega^{r}_{j}]^{{\rm{T}}}$.
3. 自适应控制器设计
对于系统(8), 在$\rho(\boldsymbol{x})$、$\varphi(\boldsymbol{x})$和$\vartheta(\boldsymbol{\chi})$未知情况下, 用前文式(18)分别对$\rho(\boldsymbol{x})$、$\varphi(\boldsymbol{x})$和$\vartheta(\boldsymbol{\chi})$进行在线逼近, 新的控制器为:
$ \begin{align} \left\{\begin{aligned} & u=\hat{\varphi}^{-1}(\boldsymbol{x})(v-\hat{\rho}(\boldsymbol{x})-\hat{\vartheta}(\boldsymbol{\chi}))\\ & v=\ddot{y}_{m}+k_{1}\dot{e}+k_{2}e \end{aligned}\right. \end{align} $
(19) 其中, $\hat{\rho}(\boldsymbol{x})$、$\hat{\varphi}(\boldsymbol{x})$和$\hat{\vartheta}(\boldsymbol{\chi})$表达式如下:
$ \begin{align} \left\{\begin{aligned} & \hat{\rho}(\boldsymbol{x}|\boldsymbol{\Theta}_{\rho})=\frac{1}{2} \boldsymbol{\Theta}_{\rho}^{{\rm{T}}}[\boldsymbol{\xi}_{\rho, l}(\boldsymbol{x})+\boldsymbol{\xi}_{\rho, r}(\boldsymbol{x})] \\ & \hat{\varphi}(\boldsymbol{x}|\boldsymbol{\Theta}_{\varphi})= \frac{1}{2}\boldsymbol{\Theta}_{\varphi}^{{\rm{T}}}[ \boldsymbol{\xi}_{\varphi, l}(\boldsymbol{x})+\boldsymbol{\xi}_{\varphi, r}(\boldsymbol{x})] \\ & \hat{\vartheta}(\boldsymbol{\chi}|\boldsymbol{\Theta}_{\vartheta})= \frac{1}{2}\boldsymbol{\Theta}_{\vartheta}^{{\rm{T}}}[\boldsymbol{\xi}_{\vartheta, l}(\boldsymbol{\chi})+\boldsymbol{\xi}_{\vartheta, r}(\chi)] \\ \end{aligned}\right. \end{align} $
(20) 定义本文二型模糊逻辑系统最优参数向量$\boldsymbol{\Theta}^{\ast}$和对应的逼近误差$\omega$分别为:
$ \begin{align} \left\{\begin{aligned} & \boldsymbol{\Theta}_{\rho}^{\ast}={\rm{arg}} \min\limits_{ \boldsymbol{\Theta}_\rho\in\Omega_{\rho}}[{\rm{sup}}_{\boldsymbol{x}\in {\bf R}^{6}}|\hat{\rho}(\boldsymbol{x}|\boldsymbol{\Theta}_{\rho})-\rho(\boldsymbol{x})|]\\ & \boldsymbol{\Theta}_{\varphi}^{\ast}={\rm{arg}} \min\limits_{ \boldsymbol{\Theta}_\varphi\in\Omega_{\varphi}}[{\rm{sup}}_{\boldsymbol{\boldsymbol{x}}\in {\bf R}^{6}}|\hat{\varphi}(\boldsymbol{x}|\boldsymbol{\Theta}_{\varphi})-\varphi(\boldsymbol{x})|]\\ & \boldsymbol{\Theta}_{\vartheta}^{\ast}={\rm{arg}} \min\limits_{ \boldsymbol{\Theta}_\vartheta\in\Omega_{\vartheta}}[{\rm{sup}}_{\boldsymbol{\chi}\in {\bf R} ^{7}}|\hat{\vartheta}(\boldsymbol{\chi}|\boldsymbol{\Theta}_{\vartheta})-\vartheta(\boldsymbol{\chi})|]\\ \end{aligned}\right. \end{align} $
(21) $ \begin{align} \omega=\, &[\hat{\rho}(\boldsymbol{x}|\boldsymbol{\Theta}_{\rho}^{\ast})-\rho(\boldsymbol{x})]+[\hat{\varphi}(\boldsymbol{x}|\boldsymbol{\Theta}_{\varphi}^{\ast})-\varphi( \boldsymbol{x})]u + \notag\\ &[\hat{\vartheta}(\boldsymbol{\chi}|\boldsymbol{\Theta}_{\vartheta}^{\ast})-\vartheta(\boldsymbol{\boldsymbol{\chi}})] \end{align} $
(22) 应用式(19)到式(8)并结合上面的定义, 经过几步直接的运算后, 可得如下的跟踪误差方程:
$ \begin{align} \dot{\boldsymbol{e}} =\, & A\boldsymbol{e}+B\Bigg[\omega+ \tilde{\boldsymbol{\Theta}}_{\varphi}^{{\rm{T}}}\frac{ \boldsymbol{\xi}_{\varphi, l}(\boldsymbol{x})+\boldsymbol{\xi}_{\varphi, r}(\boldsymbol{x})}{2} u +\notag\\ & \tilde{\boldsymbol{\Theta}}_{\rho}^{{\rm{T}}}\frac{ \boldsymbol{\xi}_{\rho, l} (\boldsymbol{x})+\boldsymbol{\xi}_{\rho, r}(\boldsymbol{x})}{2} + \tilde{\boldsymbol{\Theta}}_{\vartheta}^{{\rm{T}}} \frac{ \boldsymbol{\xi}_{\vartheta, l}(\boldsymbol{\chi})+ \boldsymbol{\xi}_{\vartheta, r}(\boldsymbol{\chi})}{2}\Bigg] \end{align} $
(23) 其中, $\tilde{\boldsymbol{\Theta}}_{\rho}=\boldsymbol{\Theta}_{\rho}-\boldsymbol{\Theta}_{\rho}^{\ast}$, $\tilde{\boldsymbol{\Theta}}_{\varphi}=\boldsymbol{\Theta}_{\varphi}-\boldsymbol{\Theta}_{\varphi}^{\ast}$, $\tilde{\boldsymbol{\Theta}}_{\vartheta}=\boldsymbol{\Theta}_{\vartheta}-\boldsymbol{\Theta}_{\vartheta}^{\ast}$, $A=\left[ {array}{cc} 0 & 1\\ -k_{2} & -k_{1} \\ {array} \right], $ $B=\left[ {array}{c} 0 \\ 1 \\ {array} \right], $ $\boldsymbol{e}=[{e}, \dot{{e}}]^{{\rm{T}}}$.
定理1. 对于系统(8)在状态变量和负载电流可测的条件下, 相应的控制器为式(19), $\boldsymbol \Theta_{\rho}$、$\boldsymbol \Theta_{\varphi}$和$\boldsymbol \Theta_{\vartheta}$的自适应律为式(24)~(26):
$ \begin{align} \dot{\boldsymbol{\Theta}}_{\rho}= \begin{cases} -\frac{\gamma_{1}}{2}\boldsymbol{e}^{{\rm{T}}}PB [ \boldsymbol{\xi}_{\rho, l}(\boldsymbol{x})+\boldsymbol{\xi}_{\rho, r}(\boldsymbol{x})], \\ \ \ \ \ \ \ \ {\rm{if}}~ ( \| \boldsymbol{\Theta}_{\rho}\|<M_{\rho}) \ {\rm{or}}~ \{\|\boldsymbol{\Theta}_{\rho}\|=M_{\rho}\\ \ \ \ \ \ \ \ {\rm{and}}\ B^{\rm{T}}P\boldsymbol{e}[ \boldsymbol{\xi}_{\rho, l}(\boldsymbol{x})+\boldsymbol{\xi}_{\rho, r}(\boldsymbol{x})]^{\rm{T}}\boldsymbol{\Theta}_{\rho} \geq0\}\\ P_{\rho}[\cdot], \ \ {\rm{if}} \ (\| \boldsymbol{\Theta}_{\rho}\|=M_{\rho}) \ {\rm{and}}\ \{B^{{\rm{T}}}P\boldsymbol{e}\times \\ \ \ \ \ \ \ \ [ \boldsymbol{\xi}_{\rho, l}(\boldsymbol{x})+\boldsymbol{\xi}_{\rho, r}(\boldsymbol{x})]^{{\rm{T}}} \boldsymbol{\Theta}_{\rho}<0\} \end{cases} \end{align} $
(24) 其中, 投影算子$P_{\rho}[\cdot]$可表示为:
$ \begin{align*} P_{\rho} [\cdot] = &-\frac{\gamma_{1}}{2}{\boldsymbol e}^{{\rm{T}}}PB \ [ \boldsymbol{\xi}_{\rho, l}(\boldsymbol{x})+\boldsymbol{\xi}_{\rho, r}(\boldsymbol{x})] +\\ & \frac{\gamma_{1}}{2}\frac{B^{{\rm{T}}}P\boldsymbol{e}[ \boldsymbol{\xi}_{\rho, l}(\boldsymbol{x})+\boldsymbol{\xi}_{\rho, r}(\boldsymbol{x})]^{{\rm{T}}}\boldsymbol{\Theta}_{\rho}}{|| \boldsymbol{\Theta}_{\rho}||^{2} }\boldsymbol{\Theta}_{\rho} \end{align*} $
$ \begin{align} \dot{\boldsymbol{\Theta}}_{\varphi}=\begin{cases} -\frac{\gamma_{2}}{2}\boldsymbol{e}^{{\rm{T}}}PB[ \boldsymbol{\xi}_{\varphi, l}(\boldsymbol{x})+\boldsymbol{\xi}_{\varphi, r}(\boldsymbol{x})]u, \\ \ \ \ \ {\rm{if}} \ (\| \boldsymbol{\Theta}_{\varphi}\|<M_{\varphi}) \ {\rm{or}}\ \{\| \boldsymbol{\Theta}_{\varphi}\|=M_{\varphi} \ \\ \ \ \ \ {\rm{and}} \ B^{{\rm{T}}}P\boldsymbol{e}[ \boldsymbol{\xi}_{\varphi, l}(\boldsymbol{x})+\boldsymbol{\xi}_{\varphi, r}(\boldsymbol{x})]^{{\rm{T}}}\boldsymbol{\Theta}_{\varphi}u \geq0\}\\ P_{\varphi}[\cdot], \ {\rm{if}}\ (\| \boldsymbol{\Theta}_{\varphi}\|=M_{\varphi}) \ {\rm{and}} \ \{B^{\rm{T}}P\boldsymbol{e}\times \\ \ \ \ \ [ \boldsymbol{\xi}_{\varphi, l}(\boldsymbol{x})+\boldsymbol{\xi}_{\varphi, r}(\boldsymbol{x})]^{\rm{T}} \boldsymbol{\Theta}_{\varphi}u<0\} \end{cases} \end{align} $
(25) 其中, 投影算子$P_{\varphi}[\cdot]$可表示为:
$ \begin{align*} P_{\varphi} [\cdot] = \, &-\frac{\gamma_{2}}{2}\boldsymbol{e}^{{\rm{T}}}PB \ [ \boldsymbol{\xi}_{\varphi, l}(\boldsymbol{x})+\boldsymbol{\xi}_{\varphi, r}(\boldsymbol{x})]u + \\ & \frac{\gamma_{2}}{2}\frac{B^{{\rm{T}}}P\boldsymbol{e}[ \boldsymbol{\xi}_{\varphi, l}(\boldsymbol{x})+\boldsymbol{\xi}_{\varphi, r}(\boldsymbol{x})]^{{\rm{T}}}\boldsymbol{\Theta}_{\varphi}u}{|| \boldsymbol{\Theta}_{\varphi}||^{2} }\boldsymbol{\Theta}_{\varphi} \ \\ \end{align*} $
$ \begin{align} \dot{\boldsymbol{\Theta}}_{\vartheta}=\begin{cases} -\frac{\gamma_{3}}{2}\boldsymbol{e}^{{\rm{T}}}PB[ \boldsymbol{\xi}_{\vartheta, l}(\boldsymbol{\chi})+\boldsymbol{\xi}_{\vartheta, r}(\boldsymbol{\chi})], \\ \ \ \ \ {\rm{if}} \ (\| \boldsymbol{\Theta}_{\vartheta}\|<M_{\vartheta}) \ {\rm{or}} \ \{\|\boldsymbol{\Theta}_{\vartheta}\|=M_{\vartheta} \ \\ \ \ \ \ {\rm{and}} \ B^{{\rm{T}}}P\boldsymbol{e}[ \boldsymbol{\xi}_{\vartheta, l}(\boldsymbol{\chi})+\boldsymbol{\xi}_{\vartheta, r}(\boldsymbol{\chi})]^{{\rm{T}}}\boldsymbol{\Theta}_{\vartheta} \geq0\}\\ P_{\vartheta}[\cdot], \ \ \ {\rm{if}} \ (\| \boldsymbol{\Theta}_{\vartheta}\|=M_{\vartheta}) \ {\rm{and}} \ \{B^{\rm{T}}P\boldsymbol{e}\times\\ \ \ \ \ [ \boldsymbol{\xi}_{\vartheta, l}(\boldsymbol{\chi})+ \boldsymbol{\xi}_{\vartheta, r}(\boldsymbol{\chi})]^{\rm{T}} \boldsymbol{\Theta}_{\vartheta}<0\} \end{cases} \end{align} $
(26) 其中, 投影算子$P_{\vartheta}[\cdot]$可表示为:
$ \begin{align*} P_{\vartheta} [\cdot] = &-\frac{\gamma_{3}}{2}\boldsymbol{e}^{\rm{T}}PB \ [ \boldsymbol{\xi}_{\vartheta, l}(\boldsymbol{\chi})+\boldsymbol{\xi}_{\vartheta, r}(\boldsymbol{\chi})]+ \\ & \frac{\gamma_{3}}{2}\frac{B^{\rm{T}}P\boldsymbol{e}[ \boldsymbol{\xi}_{\vartheta, l}(\boldsymbol{\chi})+\boldsymbol{\xi}_{\vartheta, r}(\boldsymbol{\chi})]^{\rm{T}}\boldsymbol{\Theta}_{\vartheta}}{|| \boldsymbol{\Theta}_{\vartheta}||^{2} }\boldsymbol{\Theta}_{\vartheta} \ \end{align*} $
其中, $\gamma_{1}$、$\gamma_{2}$和$\gamma_{3}$为正常数, $M_{\rho}$、$ M_{\varphi}$和$M_{\vartheta}$为系统约束集, $P=P^{T}$为满足下面的Lyapunov方程:
$ \begin{align} A^{\rm{T}}P+PA=-Q \end{align} $
(27) 那么跟踪误差$\boldsymbol{e}$是一致最终有界的.
证明. 取Lyapunov函数为:
$ \begin{align} {V}=\frac{1}{2}\boldsymbol{e}^{\rm{T}}P\boldsymbol{e}+\frac{1}{2\gamma_{1}}\tilde{\boldsymbol{\Theta}}_{\rho}^{\rm{T}}\tilde{\boldsymbol{\Theta}}_{\rho} +\frac{1}{2\gamma_{2}}\tilde{\boldsymbol{\Theta}}_{\varphi}^{\rm{T}}\tilde{\boldsymbol{\Theta}}_{\varphi} +\frac{1}{2\gamma_{3}}\tilde{\boldsymbol{\Theta}}_{\vartheta}^{\rm{T}}\tilde{\boldsymbol{\Theta}}_{\vartheta} \end{align} $
(28) 其中, ${P}={P}^{\rm T}$.微分(28)可推出:
$ \begin{align} {\dot{V}} =\, & \frac{1}{2}\dot{{\boldsymbol{e}}}^{ \rm{T}}{P}{\boldsymbol{e}}+\frac{1}{2}{\boldsymbol{e}}^{\rm{T}}{P}\dot{{\boldsymbol{e}}} + \frac{1}{\gamma_{1}}\tilde{\boldsymbol{\Theta }}_{\rho}^{\rm{T}}\dot{\tilde{\boldsymbol{\Theta }}}_{\rho} \notag+ \\ & \frac{1}{\gamma_{2}}\tilde{\boldsymbol{\Theta }}_{\varphi}^{\rm{T}}\dot{\tilde{\boldsymbol{\Theta }}}_{\varphi} + \frac{1}{\gamma_{3}}\tilde{\boldsymbol{\Theta }}_{\vartheta}^{\rm{T}}\dot{\tilde{\boldsymbol{\Theta }}}_{\vartheta} \end{align} $
(29) 基于式(23)和式(27), 可以得到:
$ \begin{align} { \dot{V} }& =-\frac{1}{2}\boldsymbol{e}^{\rm{T}}Q\boldsymbol{e} +\boldsymbol{e}^{\rm{T}}PB\Bigg[\omega+\tilde{\boldsymbol{\Theta}}_{\rho}^{\rm{T}}\frac{ \boldsymbol{\xi}_{\rho, l}(\boldsymbol{x}) +\boldsymbol{\xi}_{\rho, r}(\boldsymbol{x})}{2} + \notag\\ & \tilde{\boldsymbol{\Theta}}_{\varphi}^{\rm{T}}\frac{ \boldsymbol{\xi}_{\varphi, l}(\boldsymbol{x})+\boldsymbol{\xi}_{\varphi, r}(\boldsymbol{x})}{2} u + \tilde{\boldsymbol{\Theta}}_{\vartheta}^{\rm{T}}\frac{ \boldsymbol{\xi}_{\vartheta, l} (\boldsymbol{\chi})+\boldsymbol{\xi}_{\vartheta, r}(\boldsymbol{\chi})}{2}\Bigg] + \notag\\ & \frac{1}{\gamma_{1}}\tilde{\boldsymbol{\Theta }}_{\rho}^{\rm{T}} \dot{\tilde{\boldsymbol{\Theta }}}_{\rho} + \frac{1}{\gamma_{2}}\tilde{\boldsymbol{\Theta }}_{\varphi}^{\rm{T}} \dot{\tilde{\boldsymbol{\Theta }}}_{\varphi} + \frac{1}{\gamma_{3}}\tilde{\boldsymbol{\Theta }}_{\vartheta}^{\rm{T}} \dot{\tilde{\boldsymbol{\Theta }}}_{\vartheta} \end{align} $
(30) 上式经简单整理可得:
$ \begin{align} { \dot{V}} =\, & -\frac{1}{2}\boldsymbol{e}^{\rm{T}}Q\boldsymbol{e} +\boldsymbol{e}^{\rm{T}}PB\omega\notag+ \\ & \frac{1}{\gamma_{1}}\tilde{\boldsymbol{\Theta }}_{\rho}^{\rm{T}}\left[\dot{\tilde{\boldsymbol{\Theta }}}_{\rho}+\gamma_{1}\boldsymbol{e}^{\rm{T}}PB\frac{ \boldsymbol{\xi}_{\rho, l}(\boldsymbol{x}) +\boldsymbol{\xi}_{\rho, r}(\boldsymbol{x})}{2}\right]\notag+ \\ & \frac{1}{\gamma_{2}}\tilde{\boldsymbol{\Theta }}_{\varphi}^{\rm{T}}\left[\dot{\tilde{\boldsymbol{\Theta }}}_{\varphi}+\gamma_{2}\boldsymbol{e}^{\rm{T}}PB\frac{ \boldsymbol{\xi}_{\varphi, l}(\boldsymbol{x})+\boldsymbol{\xi}_{\varphi, r}(\boldsymbol{x})}{2}u\right]\notag+ \\ & \frac{1}{\gamma_{3}}\tilde{\boldsymbol{\Theta }}_{\vartheta}^{\rm{T}}\left[\dot{\tilde{\boldsymbol{\Theta }}}_{\vartheta}+\gamma_{3}\boldsymbol{e}^{\rm{T}}PB\frac{ \boldsymbol{\xi}_{\vartheta, l}(\boldsymbol{\chi}) +\boldsymbol{\xi}_{\vartheta, r}(\boldsymbol{\chi})}{2}\right] \end{align} $
(31) 由于$\dot{\tilde{\boldsymbol{\Theta}}}_{\rho}=\dot{\boldsymbol{\Theta}}_{\rho}$, $\dot{\tilde{\boldsymbol{\Theta}}}_{\varphi}=\dot{\boldsymbol{\Theta}}_{\varphi}$, $\dot{\tilde{\boldsymbol{\Theta}}}_{\vartheta}=\dot{\boldsymbol{\Theta}}_{\vartheta}$, 并结合式(24)~(26)可得:
$ \begin{align} {\rm{ \dot{V}}}\leq-\frac{1}{2}\boldsymbol{e}^{\rm{T}}Q\boldsymbol{e}+\boldsymbol{e}^{\rm{T}}PB\omega \end{align} $
(32) 由于矩阵$P$和矩阵$Q$为正定矩阵, 特征值均大于零, 那么以上方程满足:
$ \begin{align} {\rm{\dot{V}}} & \leq -\frac{1}{2}\lambda_{\rm{min}}(Q) \|\boldsymbol{e}\|^{2}+\|w_{0}\|\lambda_{\rm{max}}(P)\|\boldsymbol{e}\| = \notag\\ & -\frac{1}{2}\|\boldsymbol{e}\|[ \ \lambda_{\rm{min}}(Q) \|\boldsymbol{e}\|-2\lambda_{\rm{max}}(P) \|\omega_{0}\| ] \end{align} $
(33) 其中, $\omega_{0}$为逼近误差$\omega$的界, $\lambda_{{\rm{min}}}(Q)$为矩阵$Q$的最小特征值, $\lambda_{{\rm{max}}}(P) $为矩阵$P$最大特征值.
从上面方程式(33)可以看出${\dot{V}}$是变号函数, 但只要方括号中的项大于零, 那么可得${\dot{V}}<0$, 即满足下式:
$ \begin{equation} \|\boldsymbol{e}\|>\frac{2\lambda_{ \rm{max}}(P)\|\omega_{0}\|}{\lambda_{ \rm{min}}(Q)} \end{equation} $
(34) 由Lyapunov理论可知, 只要误差项大于方程式(34)的右边项, 那么跟踪误差就会减小.这也说明跟踪误差的有界性, 即:
$ \begin{equation} \|\boldsymbol{e}\|\leq\frac{2\lambda_{\rm{max}}(P)\|\omega_{0}\|}{\lambda_{\rm{min}}(Q)} \end{equation} $
(35) 由式(35)的可知, 该系统的跟踪误差大小取决于系统逼近误差的界$\omega_{0}$.由万能逼近特性可知, 自适应模糊系统的建模误差是能保证一定的精度.
4. 仿真结果与分析
把本文所提出的基于Ⅱ型模糊逻辑系统的自适应控制器用于图 1所示的PEMFC空气供给系统的过氧比控制仿真.本文设计的控制器不需要模型参数已知的条件, 只是在控制仿真验证实验时、被控对象模型需用文献[5]中的模型参数. PEMFC空气供给系统的状态空间模型为:
$ \begin{align} \left\{\begin{aligned} & \dot{\boldsymbol x}={\boldsymbol f}(\boldsymbol x)+{\boldsymbol g}\cdot u+{\boldsymbol d}\cdot {\rm I}_{st}\\ & y=h(\boldsymbol x) \end{aligned}\right. \end{align} $
(36) 在控制器设计中, 为了方便二型模糊集合统一描述和编程, 用统一向量$\boldsymbol{x}$表示.对于控制器(19)中的$\hat{\rho}(\boldsymbol{x})$和$\hat{\varphi}(\boldsymbol{x})$, 设向量$\boldsymbol{x}=[x_{1}, \cdots, x_{6}]^{\rm{T}}$; 对于$\hat{\vartheta}(\boldsymbol{\chi})$, 设向量$\boldsymbol{x}=\boldsymbol{\chi}=[x_{1}, \cdots, x_{7}]^{\rm{T}}$.
步骤1. 建立二型模糊逻辑系统
前提变量$x_{i}$经归一化处理后、取5个如图 4所示的二型模糊集合, 用$\{A_{i}^{1}, \cdots, A_{i}^{5}\}$来表示, 其相应的高斯隶属度函数选择如下: $\overline{\mu}_{\tilde{A}_{i}^{1}}=1/{\rm{exp}}[4x_{i}^{2}], \underline{\mu}_{\tilde{A}_{i}^{1}}=1/{\rm{exp}}[5x_{i}^{2}]$; $\overline{\mu}_{\tilde{A}_{i}^{2}}=1/{\rm{exp}}[4(x_{i}-0.25)^{2}], \underline{\mu}_{\tilde{A}_{i}^{2}}=1/{\rm{exp}}[5(x_{i}-0.25)^{2}]$; $\overline{\mu}_{\tilde{A}_{i}^{3}}=1/{\rm{exp}}[4(x_{i}-0.5)^{2}], \underline{\mu}_{\tilde{A}_{i}^{3}}=1/{\rm{exp}}[5(x_{i}-0.5)^{2}]$; $\overline{\mu}_{\tilde{A}_{i}^{4}}=1/{\rm{exp}}[4(x_{i}-0.75)^{2}], \underline{\mu}_{\tilde{A}_{i}^{4}}=1/{\rm{exp}}[5(x_{i}-0.75)^{2}]$; $\overline{\mu}_{\tilde{A}_{i}^{5}}=1/{\rm{exp}}[4(x_{i}-1)^{2}], \underline{\mu}_{\tilde{A}_{i}^{5}}=1/{\rm{exp}}[5(x_{i}-1)^{2}], $依据前文第2节的步骤1~步骤4可建立二型模糊逻辑系统$\hat{\rho}$、$\hat{\varphi}$和$\hat{\vartheta}$.
步骤2. 仿真参数选择
在仿真中, 被控对象的模型参数按照文献[5]的附录D和E选取; 而本文控制器的参数选择采用反复凑试法、直至出现满意的响应.状态变量的初始值基于文献[7]选取为$\boldsymbol{x}(0)=[5.4\times10^{3}, 1.5\times10^{5}, 3\times10^{-2}, 1\times10^{-3}, 8\times10^{-3}, 1.25\times10^{5}]$; 凑试自适应调整参数$\gamma_{1}=10^{2}$、$\gamma_{2}=10$和$\gamma_{3}=10^{2}$, 约束集取$M_{\rho}=10^{4}$、$M_{\varphi}=10^{2}$和$M_{\vartheta}=10^{3}$.取$k_{1}=60, k_{2}=800$, 选取正定对称矩阵$Q=\left[ {array}{cc} 16\times 10^{3} & 12\times 10^{2}\\ 12\times 10^{2} & 100 \\ {array} \right]$, 解Lyapunov方程(27)可得$ P=\left[ {array}{cc} 200 & 10\\ 10 & 1 \\ {array} \right]$为一正定对称矩阵.
步骤3. 第一种负载电流情况的仿真
为了验证基于Ⅱ型模糊系统的自适应控制有效性, 首先选取国际文献中常用图 5(a)所示的负载电流进行仿真实验.
由式(4)可知, 该电流对应的参考过氧比或者理想过氧比$\lambda_{{\rm O_{2, ref}}}$如图 5(b)实线所示. PEMFC系统过氧比$\lambda_{{\rm O_{2}}}$跟踪曲线与跟踪误差$e$曲线分别如图 5(b)虚线和图 5(c)所示, 图 5(d)为输入工作电压.由图 5(b)可知, 在图 5(a)负载电流工况下, 基于Ⅱ型模糊系统的自适应控制器能保证系统过氧比$\lambda_{{\rm O_{2}}}$很好地跟踪理想过氧比$\lambda_{{\rm O_{2}, ref}}$.由图 5(c)可知, 该控制器可保证跟踪误差一致最终有界, 从而保证了系统具有良好的跟踪性能与稳定性能.
步骤4. 第二种负载电流情况的仿真
为了进一步验证基于Ⅱ模糊系统的自适应控制器的鲁棒性, 将第一种电流切换到图 6(a)所示频繁变化的负载电流再次仿真实验. PEMFC系统过氧比$\lambda_{{\rm O_{2}}}$跟踪曲线与跟踪误差$e$曲线分别如图 6(b)和图 6(c)所示, 图 6(d)为输入工作电压.由图 6(b)可知, 在6(a)所示的变化频繁的负载电流工况下, 基于Ⅱ型模糊系统的自适应控制器同样可实现对PEMFC空气供给系统的合理控制, 使得系统过氧比$\lambda_{{\rm O_{2}}}$也能很好地跟踪理想过氧比$\lambda_{{\rm O_{2, ref}}}$.由图 6(c)可知, 在该控制器作用下, 跟踪误差也取得了一致最终有界性, 从而保证了该系统的稳定性与鲁棒性.
步骤5. 不确定参数的对比仿真
为了进一步验证基于Ⅱ模糊系统的自适应控制器的自适应特性, 将仿真模型中的原始电堆温度$T_{st}=80\, ^{\circ}{ \rm C}$切换到$T_{st}=75\, ^{\circ}{ \rm C}$.在原始温度$T_{st}=80\, ^{\circ}{ \rm C}$时, 采用精确线性化控制器(见附录C)所获得的过氧比$\lambda_{{\rm O_{2}}}$跟踪曲线和跟踪误差${ e}$曲线分别如图 7(a)和图 7(b)所示, 从图 7可知系统获得较好的跟踪效果.在控制器参数不变的情况下, 将模型中的原始温度$T_{st}=80\, ^{\circ}{ \rm C}$切换到$T_{st}=75\, ^{\circ}{ \rm C}$时, 图 8(a)和图 8(b)是精确线性化控制器作用下所获得的跟踪效果、而图 8(c)和图 8(d)是本文控制器作用下所获得的跟踪效果. 图 8(a)和(b)相比于图 7(a)和(b)说明精确线性化控制器在模型参数存在不确定的情况下、控制效果变差, 而图 8(c)和(d)说明本文所设计的控制器即使存在参数不确定, 也能获得良好的自适应性能.
图 7 在$T_{st}=80\, ^{\circ}{ \rm C}$时精确线性化控制器的仿真结果Fig. 7 The simulation results of exact linearization controller when $T_{st}=80\, ^{\circ}{ \rm C}$
图 8 在$T_{st}=75\, ^\circ{ \rm C}$时精确线性化控制器和本文所建议控制器的对比仿真结果Fig. 8 The simulation results of exact linearization controller and proposed controller when $T_{st}=75\, ^{\circ}{ \rm C}$5. 结论
本文针对PEMFC空气供给系统难以实现精准建模, 从而造成控制器的设计和系统稳定性分析比较复杂这一难题.提出基于Lyapunov稳定性的自适应控制器设计, 其中采用二型自适应模糊系统在线逼近PEMFC系统中的未建模动态, 并从Lyapunov函数中导出自适应参数, 从而保证了跟踪误差的有界性.从仿真中得出采用该控制器不但可以获得良好的跟踪性能, 而且具有良好的稳定性、鲁棒性与自适应特性.进一步的工作是在条件成熟的情况下, 将本文的控制方法应用到实际中.
附录A. 模型修正
注A1. 文献[5]中存在一些明显打印错误的地方: 1)在表 A1中, 从原公式的右边表达式显然可知注入阴极氧气流量$W_{{\rm O_{2}}, \rm in}(x_{2}, x_{3}, x_{4})$和注入阴极氮气流量$W_{{\rm N_{2}}, \rm in}(x_{2}, x_{3}, x_{4})$是变量$x_{2}, x_{4}$和$x_{5}$的函数、而非$x_{2}, x_{3}$和$x_{4}$的函数, 故公式左边需修改成$W_{{\rm O_{2}}, \rm in}(x_{2}, x_{4}, x_{5})$和$W_{{\rm N_{2}}, \rm in}(x_{2}, x_{4}, x_{5})$, 而方程右边可能是推导失误. 2)文献[5]中的状态方程(3)中的$J_{cm}$符号不存在, 由该文的附录A可知应是$J_{cp}$. 3)另外, 文献[5]状态方程(3)中的$\frac{nM_{{\rm O_{2}}}}{4F}$应是$-\frac{nM_{{\rm O_{2}}}}{4F}$. 4)在文献[5]中的$B_{43}-B_{45}$没有定义且该文也没有用到, 但上表中的流出阴极空气流量却用到$B_{46}-B_{49}$且在该文附录F中也没有此定义、应是$B_{1}-B_{42}$当中的某些量, 显然是误写.
表 A1 原公式和修正后公式的对比Table A1 Comparison of original formulas and revised formulas物理意义 原公式 修正后公式 注入阴极氧气流量 ${array}{c}W_{{\rm O_{2}}, \rm in}(x_{2}, x_{3}, x_{4})=((x_{2}-B_{32}-B_{33}-\\x_{5}B_{34}-x_{4}B_{35})\times(x_{2}-x_{2}B_{6})^{-1}+ \\(x_{2}B_{36}-B_{37}-x_{5}B_{38}-\\x_{4}B_{39}))e(x_{2})k(x_{2}){array}$ ${array}{c}W_{{\rm O_{2}}, \rm in}(x_{2}, x_{4}, x_{5})=((x_{2}B_{32}-B_{33}-\\x_{5}B_{34}-x_{4}B_{35})\times(x_{2}-x_{2}B_{6})^{-1}+\\(x_{2}B_{36}-K_{sm, \rm out}B_{37}-x_{5}K_{sm, \rm out}B_{38}-\\x_{4}K_{sm, \rm out}B_{39}))e(x_{2})k(x_{2}){array}$ 流出阴极空气流量 ${array}{c}W_{ca, \rm out}(x_{4}, x_{5}, x_{6})=B_{47}+x_{5}B_{48}+\\x_{4}B_{49}-x_{6}B_{46} {array}$ ${array}{c}W_{ca, \rm out}(x_{4}, x_{5}, x_{6})=B_{20}+x_{5}B_{21}+\\ x_{4}B_{22}-x_{6}B_{19}{array}$ 注入阴极氮气流量 ${array}{c}W_{{\rm N_{2}}, \rm in}(x_{2}, x_{3}, x_{4})=((x_{2}B_{23}-B_{24}-\\x_{5}B_{25}-x_{4}B_{26})\times(x_{2}-x_{2}B_{6})^{-1}+\\(x_{2}B_{27}-B_{28}-x_{5}B_{29}-\\x_{4}B_{30}))e(x_{2})k(x_{2}){array}$ ${array}{c}W_{{\rm N_{2}}, \rm in}(x_{2}, x_{4}, x_{5})=((x_{2}B_{23}-B_{24}-\\x_{5}B_{25}-x_{4}B_{26})\times(x_{2}-x_{2}B_{6})^{-1}+\\(x_{2}B_{27}-K_{sm, \rm out}B_{28}-x_{5}K_{sm, \rm out}B_{29}-\\x_{4}K_{sm, \rm out}B_{30}))e(x_{2})k(x_{2}){array}$ 流出阴极氧气流量 ${array}{c} W_{{\rm O_{2}}, \rm out}(x_{4}, x_{5}, x_{6})=-x_{4}(B_{10}-x_{5}B_{11} +\\ x_{4}B_{12}-x_{6}B_{9})\times j(x_{4}, x_{5})x_{4 }^{-1} \times\\(j(x_{4}, x_{5})B_{40}-M_{N_{2}})^{-1} \times m(x_{4}, x_{5}){array}$ ${array}{c}W_{{\rm O_{2}}, \rm out}(x_{4}, x_{5}, x_{6})=x_{4}(B_{10}+x_{5}B_{11} +\\x_{4}B_{12}-x_{6}B_{9}) \times j(x_{4}, x_{5})x_{4 }^{-1} \times\\(j(x_{4}, x_{5})B_{40}+M_{N_{2}})^{-1} \times m(x_{4}, x_{5}){array}$ 空压机驱动力矩 ${array}{c}\tau_{cm}(u, x_{1})=\frac{\eta_{cm}K_{t}(u-K_{v}x_{1})}{R_{cm}J_{cp}}{array}$ ${array}{c}\tau_{cm}(u, x_{1})=\frac{\eta_{cm}K_{t}(u-K_{v}x_{1})}{R_{cm}}{array} $ 空压机负载力矩 ${array}{c}\tau_{cp}(x_{1}, x_{2})=\frac{C_{p}T_{atm}n(x_{2})W_{cp}(x_{1}, x_{2})}{\eta_{cp}J_{cp}x_{1}}{array}$ ${array}{c}\tau_{cp}(x_{1}, x_{2})=\frac{C_{p}T_{atm}n(x_{2})W_{cp}(x_{1}, x_{2})}{\eta_{cp}x_{1}}{array}$ 下面对以上面第一个表达式注入阴极氧气流量$W_{{\rm O_{2}}, \rm in}(x_{2}, x_{4}, x_{5})$为例说明其修正过程:
由文献[7]可知, $W_{{\rm O_{2}}, \rm in}(x_{2}, x_{4}, x_{5})$可表示为:
$ \begin{align} \begin{cases} W_{{\rm O_{2}}, \rm in}(x_{2}, x_{4}, x_{5})=X_{{\rm O_{2}}, ca, \rm in}W_{a, ca, in}(x_{2}, x_{4}, x_{5})\\ W_{a, ca, in}(x_{2}, x_{4}, x_{5})=\frac{W_{ca, \rm in}(x_{2}, x_{4}, x_{5})}{1+\omega_{ca, \rm in}(x_{2})} \end{cases} \end{align} $
(A1) 其中, $X_{{\rm O_{2}}, ca, \rm in}$为注入阴极空气中的氧气质量分数, $W_{a, ca, \rm in}(x_{2}, x_{4}, x_{5})$为注入阴极干燥空气流量, $\omega_{ca, \rm in}(x_{2})$为注入阴极空气的湿度, $W_{ca, \rm in}(x_{2}, x_{4}, x_{5})$为注入阴极空气流量.首先, $X_{{\rm O_{2}}, ca, \rm in}$和$\omega_{ca, \rm in}(x_{2})$可分别表示为:
$ \begin{align} \begin{cases} X_{{\rm O_{2}}, ca, \rm in}=\frac{y_{{\rm O_{2}}, ca, \rm in}M_{{\rm O_{2}}}}{y_{{\rm O_{2}}, ca, \rm in}M_{{\rm O_{2}}}+(1-y_{{\rm O_{2}}, ca, \rm in})M_{{\rm N_{2}}}}\\ \omega_{ca, in}(x_{2})=\frac{M_{v}}{y_{{\rm O_{2}}, ca, \rm in}M_{{\rm O_{2}}}+(1-y_{{\rm O_{2}}, ca, \rm in}) M_{{\rm N_{2}}}}\times \\ \qquad \qquad \qquad \frac{p_{v, ca, \rm in}}{p_{a, ca, \rm in}(x_{2})} \end{cases} \end{align} $
(A2) 其中, $M_{{\rm O_{2}}}$、$M_{{\rm N_{2}}}$和$M_{v}$分别为氧气、氮气和水蒸气的摩尔质量, $y_{{\rm O_{2}}, ca, \rm in}$为注入阴极空气中氧气的摩尔质量分数, $p_{v, ca, \rm in}$和$p_{a, ca, \rm in}(x_{2})$分别为注入阴极的水蒸气压强和干燥空气的压强.接下来, $p_{v, ca, \rm in}$和$p_{a, ca, \rm in}(x_{2})$可分别表示为:
$ \begin{align} \begin{cases} p_{v, ca, \rm in}=\phi_{ca, \rm in}p_{sat, T_{cl}}\\ p_{a, ca, \rm in}(x_{2})=x_{2}-p_{v, sm}(x_{2})+p_{\Delta, hm} \end{cases} \end{align} $
(A3) 其中, $\phi_{ca, \rm in}$为注入阴极空气的相对湿度, $p_{sat, T_{cl}}$为经冷却处理后水蒸气的饱和压强, $p_{v, sm}(x_{2})$为供给管道内水蒸气压强, $p_{\Delta, hm}$为加湿处理后的水蒸气压强与注入阴极的水蒸气压强之差.再有, $p_{v, sm}(x_{2})$和$p_{\Delta, hm}$可分别表示为:
$ \begin{align} \begin{cases} p_{v, sm}(x_{2})=x_{2}\frac{\phi_{atm}p_{sat, T_{atm}}}{p_{atm}}\notag\\ p_{\Delta, hm}=\phi_{des}p_{sat, T_{cl}}-\phi_{ca, \rm in}p_{sat, T_{cl}} \end{cases} \end{align} $
其中, $p_{atm}$、$\phi_{atm}$和$p_{sat, T_{atm}}$分别为标准状况下空气压强、空气的相对湿度和空气中水蒸气饱和压强, $\phi_{des}$为加湿处理后的空气相对湿度.
由文献[7] (在文献[5]中的参考文献[35])可知, 注入阴极空气流量$W_{ca, \rm in}(x_{2}, x_{4}, x_{5})$可表示为:
$ \begin{align} W_{ca, \rm in}(x_{2}, x_{4}, x_{5})=\, &W_{a, hm}(x_{2}, x_{4}, x_{5})+ \notag\\ &W_{v, hm}(x_{2}, x_{4}, x_{5}) \end{align} $
(A4) 其中, $W_{a, hm}(x_{2}, x_{4}, x_{5})$和$W_{v, hm}(x_{2}, x_{4}, x_{5})$分别为加湿处理后干燥空气流量与水蒸气流量, 可分别表示为:
$ \begin{align} \begin{cases} W_{a, hm}(x_{2}, x_{4}, x_{5})=\frac{W_{sm, \rm out}(x_{2}, x_{4}, x_{5})}{1+\omega_{sm}(x_{2})}\\ W_{v, hm}(x_{2}, x_{4}, x_{5})=\frac{M_{v}}{M_{a}}\frac{\phi_{des}p_{sat, T_{cl}}}{x_{2}-p_{v, sm}(x_{2})}\times \\ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ W_{a, hm}(x_{2}, x_{4}, x_{5}) \end{cases} \end{align} $
(A5) 其中, $\omega_{sm}(x_{2})$为供给管道内空气湿度, $W_{sm, \rm out}(x_{2}, x_{4}, x_{5})$为流出供给管道的空气流量, 可分别表示为:
$ \begin{align} \begin{cases} \omega_{sm}(x_{2})=\frac{M_{v}}{M_{a}}\frac{p_{v, sm}(x_{2})}{x_{2}-p_{v, sm}(x_{2})}\\ W_{sm, \rm out}(x_{2}, x_{4}, x_{5})=K_{sm, \rm out}(x_{2}-p_{ca}(x_{4}, x_{5})) \end{cases} \end{align} $
(A6) 其中, $K_{sm, \rm out}$为供给管道孔口常数, $p_{ca}(x_{4}, x_{5})$为阴极压强. $p_{ca}(x_{4}, x_{5})$可表示为:
$ \begin{align} \begin{cases} p_{ca}(x_{4}, x_{5})=p_{{\rm O_{2}}}(x_{4})+p_{{\rm N_{2}}}(x_{5})+p_{v, ca}\\ p_{O_{2}, ca}(x_{4})=\frac{x_{4}R_{{\rm O_{2}}}T_{st}}{V_{ca}}\\ p_{{\rm N_{2}}, ca}(x_{5})=\frac{x_{5}R_{{\rm N_{2}}}T_{st}}{V_{ca}} \end{cases} \end{align} $
(A7) 其中, $p_{v, ca}$、$p_{{\rm O_{2}}, ca}(x_{4})$和$p_{{\rm N_{2}}, ca}(x_{5})$分别为阴极内水蒸气、氧气和氮气压强, $R_{{\rm O_{2}}}$和$R_{{\rm N_{2}}}$分别为氧气和氮气的气体常数, $V_{ca}$为阴极体积.
运用式(A2)~(A7)到式(A1), 经简单整理后可得注入阴极氧气的流量$W_{{\rm O_{2}}, \rm in}(x_{2}, x_{4}, x_{5})$为:
$ \begin{align} &W_{{\rm O_{2}}, \rm in}(x_{2}, x_{4}, x_{5}) =K_{sm, \rm out}\Big(x_{2}-p_{v, ca}-\notag\\ &\qquad \frac{x_{4}R_{{\rm O_{2}}}T_{st}} {V_{ca}}- \frac{x_{5}R_{{\rm N_{2}}}T_{st}}{V_{ca}}\Big)\times\notag\\ &\qquad \frac{y_{{\rm O_{2}}, ca, \rm in}M_{{\rm O_{2}}}}{y_{{\rm O_{2}}, ca, \rm in}M_{{\rm O_{2}}}+(1-y_{{\rm O_{2}}, ca, \rm in})M_{{\rm N_{2}}}}\times \notag\\ &\qquad\left(1+\frac{M_{v}}{M_{a}}\frac{\phi_{des}p_{sat, T_{cl}}}{x_{2}- x_{2}\phi_{atm}p_{sat, T_{atm}}p_{atm}^{-1}}\right) \times \notag\\ &\qquad\left(1+\frac{M_{v}}{M_{a}}\frac{\phi_{atm}p_{sat, T_{atm}}} {x_{2}-x_{2}\phi_{atm}p_{sat, T_{atm}}p_{atm}^{-1}}\right)^{-1}\times \notag\\ &\qquad\Big(1+\frac{M_{v}}{y_{{\rm O_{2}}, ca, \rm in}M_{{\rm O_{2}}}+ (1-y_{{\rm O_{2}}, ca, \rm in})M_{{\rm N_{2}}}}\times \notag\\ & \qquad\frac{\phi_{ca, \rm in}p_{sat, T_{cl}}}{x_{2}-p_{v, sm}(x_{2})+\phi_{des}p_{sat, T_{cl}}-\phi_{ca, \rm in}p_{sat, T_{cl}}}\Big)^{-1} \end{align} $
(A8) 结合文献[5]附录C中$e(x_{2})$、$k(x_{2})$函数的定义和附录F常量$B_{32}\sim B_{39}$的定义, 式(A8)可表示为:
$ \begin{align} &W_{{\rm O_{2}}, \rm in}(x_{2}, x_{4}, x_{5})=((x_{2}B_{32}-B_{33}-x_{5}B_{34}- \notag\\ &\quad x_{4}B_{35})(x_{2}-x_{2}B_{6})^{-1}+(x_{2}B_{36}-\notag\\ &\quad K_{sm, \rm out}B_{37} - x_{5}K_{sm, \rm out}B_{38}--\notag\\ &\quad x_{4}K_{sm, \rm out}B_{39}))e(x_{2})k(x_{2}) \end{align} $
(A9) 由于篇幅原因, 其余公式更正过程不再赘述.
附录B. 状态空间函数
$ \begin{align*} &f_{1}(x_{1}, x_{2})=\frac{1}{J_{cp}}[\tau_{cm}(u, x_{1})- \\&\qquad\tau_{cp}(x_{1}, x_{2})-\frac{\eta_{cm}K_{t}} {R_{cm}}u] \\ &f_{2}(x_{1}, x_{2}, x_{3}, x_{4}, x_{5})=\frac{\gamma R_{a}}{V_{sm}}[T_{cp}(x_{2})W_{cp}(x_{1}, x_{2}) - \notag\\ &\qquad T_{sm}(x_{2}, x_{3})W_{sm, {\rm out}}(x_{2}, x_{4}, x_{5})]\notag\\ &f_{3}(x_{1}, x_{2}, x_{4}, x_{5})=W_{cp}(x_{1}, x_{2}) - \\ &\qquad W_{sm, {\rm out}}(x_{2}, x_{4}, x_{5})\notag\\ &f_{4}(x_{2}, x_{4}, x_{5}, x_{6})=W_{{\rm O_{2}}, \rm in}(x_{2}, x_{4}, x_{5})\notag- \\ &\qquad W_{{\rm O_{2}}, \rm out}(x_{4}, x_{5}, x_{6}) \\ \, &f_{5}(x_{2}, x_{4}, x_{5}, x_{6})=W_{{\rm N_{2}}, \rm in}(x_{2}, x_{4}, x_{5})\notag- \\ &\qquad W_{{\rm N_{2}}, \rm out}(x_{4}, x_{5}, x_{6}) \notag\\ \, &f_{6}(x_{4}, x_{5}, x_{6}) =\frac{R_{a}T_{fc}} {V_{rm}M_{a}}[ W_{ca, \rm out}(x_{4}, x_{5}, x_{6})\notag- \\ &\qquad W_{rm, \rm out}(x_{6})] \end{align*} $
其中, 在参考文献[5]中的函数$T_{cp}(x_{2})$, $W_{cp}$$(x_{1}, x_{2})$, $T_{sm}$ $(x_{2}, x_{3})$, $W_{sm, \rm out}$ $(x_{2}, x_{4}, x_{5})$, $W_{{\rm N_{2}}, \rm out}$\\ $(x_{2}, x_{4}, x_{5})$和$W_{rm, \rm out}(x_{6})$的推导与整理均正确无误.但函数$W_{{\rm O_{2}}, \rm in}(x_{2}, x_{4}, x_{5})$, $W_{{\rm O_{2}}, \rm out}(x_{4}, x_{5}, x_{6})$, \\ $W_{{\rm N_{2}}, \rm in}(x_{2}, x_{4}, x_{5})$, $\tau_{cm}(u, x_{1})$, $\tau_{cp}(x_{1}, x_{2})$和$W_{ca, \rm out}(x_{4}, x_{5}, x_{6})$的推导与整理存在错误.本文结合文献[7]中的基本物理方程对以上错误函数进行复合修正与整理, 修正结果见前面的附录A.
附录C. 精确线性化控制器求解
仿真对比实验需要用到精确线性化控制器.由于系统模型的复杂性、很难手工求解, 为此需要借助软件计算工具求解.精确线性化控制器中的$\rho(\boldsymbol x)$、$\varphi(\boldsymbol x)$和$\vartheta(\boldsymbol \chi)$函数可用Matlab软件中的mupad工具箱求解, 其计算过程指令如下:
步骤1. 在Matlab命令窗口输入“mupad”指令、打开mupad工具箱, 结合本文与文献[5]附录将正确的${\boldsymbol f}(\boldsymbol x)$、${\boldsymbol g}$、${\boldsymbol d}$和$h(\boldsymbol x)$函数编写到mupad工作空间.
步骤2. 首先通过以下微分指令函数, 求解式(6)中$L_{{\boldsymbol f}}h(\boldsymbol x)$和$L_{{\boldsymbol d}}h(\boldsymbol x)$, 即:
$ \begin{align*} \left\{\begin{aligned} & L_{{\boldsymbol f}}h(\boldsymbol x)={\rm diff}(h(\boldsymbol x), x_{2})f_{2}(\boldsymbol x)+{ \rm diff}(h(\boldsymbol x), x_{4})\times \notag\\ &\ \ \ \ \ \ \ \ \ \ f_{4}(\boldsymbol x)+{\rm diff}(h(\boldsymbol x), x_{5})f_{5}(\boldsymbol x) \notag\\ & L_{{\boldsymbol d}}h(\boldsymbol x)={\rm diff}(h(\boldsymbol x), x_{4})d_{4}\notag\\ \end{aligned}\right. \end{align*} $
步骤3. 接下来通过以下微分指令函数, 求解式(7)中$L_{{\boldsymbol f}}^{2}h(\boldsymbol x)$、$L_{{\boldsymbol g}}L_{{\boldsymbol f}}h(\boldsymbol x)$和$L_{{\boldsymbol d}}L_{{\boldsymbol f}}h(\boldsymbol x)$, 即:
$ \begin{align*} \left\{\begin{aligned} &L_{{\boldsymbol f}}^{2}h(\boldsymbol x)=\sum\limits_{i=1}^{\rm 6}{\rm diff}(L_{{\boldsymbol f}}h(\boldsymbol x), x_{i})f_{i}(\boldsymbol x)\\ &L_{{\boldsymbol g}}L_{{\boldsymbol f}}h(\boldsymbol x)={\rm diff}(L_{{\boldsymbol f}}h(\boldsymbol x), x_{1})g_{1}\\ &L_{{\boldsymbol d}}L_{{\boldsymbol f}}h(\boldsymbol x)={\rm diff}(L_{{\boldsymbol f}}h(\boldsymbol x), x_{4})d_{4} \end{aligned}\right. \end{align*} $
其中, $g_{1}$和$d_{4}$为向量$g$和$d$的分量.
步骤4. 文献[11]认为实际中的${\rm I}_{st}$电流是缓慢变化的、可视作常数处理, 用以下微分指令求解$L_{{\boldsymbol d}}h(\boldsymbol x)$关于时间$t$的导数, 即:
$ \begin{align*} \frac{{\rm d}}{{\rm d}t}(L_{{\boldsymbol d}}h(\boldsymbol x))&=\sum\limits_{i=1}^{\rm 6}{\rm diff}(L_{{\boldsymbol d}}h(\boldsymbol x), x_{i})\dot{x}_{i}(\boldsymbol x) \end{align*} $
步骤5. 最后综合步骤1到步骤4, 精确线性化控制器中的$\rho(\boldsymbol x)$、$\varphi(\boldsymbol x)$和$\vartheta(\boldsymbol \chi)$函数表达式可表示如下:
$ \begin{align*} \left\{\begin{aligned} &\rho(\boldsymbol x)=L_{{\boldsymbol f}}^{2}h(\boldsymbol x)\\ &\varphi(\boldsymbol x)=L_{{\boldsymbol g}}L_{{\boldsymbol f}}h(\boldsymbol x)\\ &\vartheta(\boldsymbol \chi)=L_{{\boldsymbol d}}L_{{\boldsymbol f}}h(\boldsymbol x){\rm I}_{st}+\frac{\rm {d}}{{\rm d}t}(L_{{\boldsymbol d}}h(\boldsymbol x)){\rm I}_{st} \end{aligned}\right. \end{align*} $
-
图 2 PEMFC系统负载电流、过氧比和输出净功率的关系
Fig. 2 Power relationship of PEMFC system load current, OER, and output net
图 7 在$T_{st}=80\, ^{\circ}{ \rm C}$时精确线性化控制器的仿真结果
Fig. 7 The simulation results of exact linearization controller when $T_{st}=80\, ^{\circ}{ \rm C}$
图 8 在$T_{st}=75\, ^\circ{ \rm C}$时精确线性化控制器和本文所建议控制器的对比仿真结果
Fig. 8 The simulation results of exact linearization controller and proposed controller when $T_{st}=75\, ^{\circ}{ \rm C}$
表 1 PEMFC空气供给系统状态变量
Table 1 State variables of PEMFC air supply system
状态变量 符号 单位 空压机转速 $x_{1}=\omega_{cp}$ rad/s 供给管道内空气压强 $x_{2}=P_{sm}$ Pa 供给管道内空气质量 $x_{3}=m_{sm}$ kg 阴极内氧气质量 $x_{4}=m_{{\rm O_{2}}}$ kg 阴极内氮气质量 $x_{5}=m_{{\rm N_{2}}}$ kg 回流管道内空气压强 $x_{6}=P_{rm}$ Pa 
下载: 导出CSV
A1 原公式和修正后公式的对比
A1 Comparison of original formulas and revised formulas
物理意义 原公式 修正后公式 注入阴极氧气流量 $\begin{array}{c}W_{{\rm O_{2}}, \rm in}(x_{2}, x_{3}, x_{4})=((x_{2}-B_{32}-B_{33}-\\x_{5}B_{34}-x_{4}B_{35})\times(x_{2}-x_{2}B_{6})^{-1}+ \\(x_{2}B_{36}-B_{37}-x_{5}B_{38}-\\x_{4}B_{39}))e(x_{2})k(x_{2})\end{array}$ $\begin{array}{c}W_{{\rm O_{2}}, \rm in}(x_{2}, x_{4}, x_{5})=((x_{2}B_{32}-B_{33}-\\x_{5}B_{34}-x_{4}B_{35})\times(x_{2}-x_{2}B_{6})^{-1}+\\(x_{2}B_{36}-K_{sm, \rm out}B_{37}-x_{5}K_{sm, \rm out}B_{38}-\\x_{4}K_{sm, \rm out}B_{39}))e(x_{2})k(x_{2})\end{array}$ 流出阴极空气流量 $\begin{array}{c}W_{ca, \rm out}(x_{4}, x_{5}, x_{6})=B_{47}+x_{5}B_{48}+\\x_{4}B_{49}-x_{6}B_{46} \end{array}$ $\begin{array}{c}W_{ca, \rm out}(x_{4}, x_{5}, x_{6})=B_{20}+x_{5}B_{21}+\\ x_{4}B_{22}-x_{6}B_{19}\end{array}$ 注入阴极氮气流量 $\begin{array}{c}W_{{\rm N_{2}}, \rm in}(x_{2}, x_{3}, x_{4})=((x_{2}B_{23}-B_{24}-\\x_{5}B_{25}-x_{4}B_{26})\times(x_{2}-x_{2}B_{6})^{-1}+\\(x_{2}B_{27}-B_{28}-x_{5}B_{29}-\\x_{4}B_{30}))e(x_{2})k(x_{2})\end{array}$ $\begin{array}{c}W_{{\rm N_{2}}, \rm in}(x_{2}, x_{4}, x_{5})=((x_{2}B_{23}-B_{24}-\\x_{5}B_{25}-x_{4}B_{26})\times(x_{2}-x_{2}B_{6})^{-1}+\\(x_{2}B_{27}-K_{sm, \rm out}B_{28}-x_{5}K_{sm, \rm out}B_{29}-\\x_{4}K_{sm, \rm out}B_{30}))e(x_{2})k(x_{2})\end{array}$ 流出阴极氧气流量 $\begin{array}{c} W_{{\rm O_{2}}, \rm out}(x_{4}, x_{5}, x_{6})=-x_{4}(B_{10}-x_{5}B_{11} +\\ x_{4}B_{12}-x_{6}B_{9})\times j(x_{4}, x_{5})x_{4 }^{-1} \times\\(j(x_{4}, x_{5})B_{40}-M_{N_{2}})^{-1} \times m(x_{4}, x_{5})\end{array}$ $\begin{array}{c}W_{{\rm O_{2}}, \rm out}(x_{4}, x_{5}, x_{6})=x_{4}(B_{10}+x_{5}B_{11} +\\x_{4}B_{12}-x_{6}B_{9}) \times j(x_{4}, x_{5})x_{4 }^{-1} \times\\(j(x_{4}, x_{5})B_{40}+M_{N_{2}})^{-1} \times m(x_{4}, x_{5})\end{array}$ 空压机驱动力矩 $\begin{array}{c}\tau_{cm}(u, x_{1})=\frac{\eta_{cm}K_{t}(u-K_{v}x_{1})}{R_{cm}J_{cp}}\end{array}$ $\begin{array}{c}\tau_{cm}(u, x_{1})=\frac{\eta_{cm}K_{t}(u-K_{v}x_{1})}{R_{cm}}\end{array} $ 空压机负载力矩 $\begin{array}{c}\tau_{cp}(x_{1}, x_{2})=\frac{C_{p}T_{atm}n(x_{2})W_{cp}(x_{1}, x_{2})}{\eta_{cp}J_{cp}x_{1}}\end{array}$ $\begin{array}{c}\tau_{cp}(x_{1}, x_{2})=\frac{C_{p}T_{atm}n(x_{2})W_{cp}(x_{1}, x_{2})}{\eta_{cp}x_{1}}\end{array}$
下载: 导出CSV
-
[1] Gruber J K, Bordons C, Oliva A. Nonlinear MPC for the airflow in a PEM fuel cell using a volterra series model. Control Engineering Practice, 2012, 20(2):205-217 doi: 10.1016/j.conengprac.2011.10.014 [2] Ramos-Paja C A, Giral R, Martinez-Salamero L, Romano J, Romero A, Spagnuolo G. A PEM fuel-cell model featuring oxygen-excess-ratio estimation and power-electronics interaction. IEEE Transactions on Industrial Electronics, 2010, 57(6):1914-1924 doi: 10.1109/TIE.2009.2026363 [3] Hayati M R, Khayatian A, Dehghani M. Simultaneous optimization of net power and enhancement of PEM fuel cell lifespan using extremum seeking and sliding mode control techniques. IEEE Transactions on Energy Conversion, 2016, 31(2):688-696 http://www.wanfangdata.com.cn/details/detail.do?_type=perio&id=c0aa01670a3a8ef4e16917ce0adf032a [4] Jayakumar A, Chalmers A, Lie T T. Review of prospects for adoption of fuel cell electric vehicles in New Zealand. IET Electrical Systems in Transportation, 2017, 7(4):259-266 doi: 10.1049/iet-est.2016.0078 [5] Pilloni A, Pisano A, Usai E. Observer-based air excess ratio control of a PEM fuel cell system via high-order sliding mode. IEEE Transactions on Industrial Electronics, 2015, 62(8):5236-5246 doi: 10.1109/TIE.2015.2412520 [6] Zhou D M, Gao F, Breaz E, Ravey A, Miraoui A, Zhang K. Dynamic phenomena coupling analysis and modeling of proton exchange membrane fuel cells. IEEE Transactions on Energy Conversion, 2016, 31(4):1399-1412 doi: 10.1109/TEC.2016.2587162 [7] Pukrushpan J T, Stefanopoulou A G, Peng H. Control of Fuel Cell Power Systems:Principles, Modeling, Analysis and Feedback Design. London:Springer-Verlag, 2004. [8] Matraji I, Laghrouche S, Jemei S, Wack M. Robust control of the PEM fuel cell air-feed system via sub-optimal second order sliding mode. Applied Energy, 2013, 104:945-957 doi: 10.1016/j.apenergy.2012.12.012 [9] Ki Na W, Gou B. Feedback-linearization-based nonlinear control for PEM fuel cells. IEEE Transactions on Energy Conversion, 2008, 23(1):179-190 http://www.wanfangdata.com.cn/details/detail.do?_type=perio&id=f9367cef535a36276ef8eb9f13f31baf [10] Laghrouche S, Liu J X, Ahmed F S, Harmouche M, Wack M. Adaptive second-order sliding mode observer-based fault reconstruction for PEM fuel cell air-feed system. IEEE Transactions on Control Systems Technology, 2015, 23(3):1098-1109 doi: 10.1109/TCST.2014.2361869 [11] Talj R J, Hissel D, Ortega R, Becherif M, Hilairet M. Experimental validation of a PEM fuel-cell reduced-order model and a moto-compressor higher order sliding-mode control. IEEE Transactions on Industrial Electronics, 2010, 57(6):1906-1913 doi: 10.1109/TIE.2009.2029588 [12] Xu L F, Hu J M, Cheng S L, Fang C, Li J Q, Ouyang M G, et al. Robust control of internal states in a polymer electrolyte membrane fuel cell air-feed system by considering actuator properties. International Journal of Hydrogen Energy, 2017, 42(18):13171-13191 doi: 10.1016/j.ijhydene.2017.03.191 [13] 莫红, 王飞跃, 肖志权, 陈茜.基于区间二型模糊集合的语言动力系统稳定性.自动化学报, 2011, 37(8):1018-1024 http://www.aas.net.cn/CN/abstract/abstract17522.shtmlMo Hong, Wang Fei-Yue, Xiao Zhi-Quan, Chen Qian. Stabilities of linguistic dynamic systems based on interval type-2 fuzzy sets. Acta Automatica Sinica, 2011, 37(8):1018-1024 http://www.aas.net.cn/CN/abstract/abstract17522.shtml [14] Saha A, Konar A, Nagar A K. EEG analysis for cognitive failure detection in driving using type-2 fuzzy classifiers. IEEE Transactions on Emerging Topics in Computational Intelligence, 2017, 1(6):437-453 doi: 10.1109/TETCI.2017.2750761 [15] Sun D, Liao Q F, Ren H L. Type-2 fuzzy modeling and control for bilateral teleoperation system with dynamic uncertainties and time-varying delays. IEEE Transactions on Industrial Electronics, 2018, 65(1):447-459 doi: 10.1109/TIE.2017.2719604 [16] Sarabakha A, Fu C H, Kayacan E, Kumbasar T. Type-2 fuzzy logic controllers made even simpler:from design to deployment for UAVs. IEEE Transactions on Industrial Electronics, 2018, 65(6):5069-5077 doi: 10.1109/TIE.2017.2767546 [17] Andreu-Perez J, Cao F, Hagras H, Yang G Z. A self-adaptive online brain-machine interface of a humanoid robot through a general type-2 fuzzy inference system. IEEE Transactions on Fuzzy Systems, 2018, 26(1):101-116 doi: 10.1109/TFUZZ.2016.2637403 [18] Wang L X. A new look at type-2 fuzzy sets and type-2 fuzzy logic systems. IEEE Transactions on Fuzzy Systems, 2017, 25(3):693-706 doi: 10.1109/TFUZZ.2016.2543746 [19] 王飞跃, 莫红.关于二型模糊集合的一些基本问题.自动化学报, 2017, 43(7):1114-1141 http://www.aas.net.cn/CN/abstract/abstract19087.shtmlWang Fei-Yue, Mo Hong. Some fundamental issues on type-2 fuzzy sets. Acta Automatica Sinica, 2017, 43(7):1114-1141 http://www.aas.net.cn/CN/abstract/abstract19087.shtml [20] Wang L X. Stable adaptive fuzzy control of nonlinear systems. IEEE Transactions on Fuzzy Systems, 1993, 1(2):146-155 doi: 10.1109/91.227383 [21] 王永富, 王殿辉, 柴天佑.一个具有完备性和鲁棒性的模糊规则提取算法.自动化学报, 2010, 36(9):1337-1342 http://www.aas.net.cn/CN/abstract/abstract17330.shtmlWang Yong-Fu, Wang Dian-Hui, Chai Tian-You. Extraction of fuzzy rules with completeness and robustness. Acta Automatica Sinica, 2010, 36(9):1337-1342 http://www.aas.net.cn/CN/abstract/abstract17330.shtml [22] Wang Y F, Wang D H, Chai T Y. Extraction and adaptation of fuzzy rules for friction modeling and control compensation. IEEE Transactions on Fuzzy Systems, 2011, 19(4):682-693 doi: 10.1109/TFUZZ.2011.2134104 [23] Wu T S, Karkoub M, Wang H W, Chen H S, Chen T H. Robust tracking control of MIMO underactuated nonlinear systems with dead-zone band and delayed uncertainty using an adaptive fuzzy control. IEEE Transactions on Fuzzy Systems, 2017, 25(4):905-918 doi: 10.1109/TFUZZ.2016.2586970 [24] Cervantes J, Yu W, Salazar S, Chairez I. Takagi-Sugeno dynamic neuro-fuzzy controller of uncertain nonlinear systems. IEEE Transactions on Fuzzy Systems, 2017, 25(6):1601-1615 doi: 10.1109/TFUZZ.2016.2612697 [25] Zhang X Y, Xu Z S, Su C Y, Li Z, Li X M, Xiong C H, et al. Fuzzy approximator based adaptive dynamic surface control for unknown time delay nonlinear systems with input asymmetric hysteresis nonlinearities. IEEE Transactions on Systems, Man, and Cybernetics:Systems, 2017, 47(8):2218-2232 doi: 10.1109/TSMC.2016.2641926 [26] Mendel J M. On KM algorithms for solving type-2 fuzzy set problems. IEEE Transactions on Fuzzy Systems, 2013, 21(3):426-446 doi: 10.1109/TFUZZ.2012.2227488 [27] Liang Q L, Mendel J M. MPEG VBR video traffic modeling and classification using fuzzy technique. IEEE Transactions on Fuzzy Systems, 2001, 9(1):183-193 doi: 10.1109/91.917124 [28] Mendel J M, John R I B. Type-2 fuzzy sets made simple. IEEE Transactions on Fuzzy Systems, 2002, 10(2):117-127 doi: 10.1109/91.995115 [29] Mendel J M, John R I, Liu F L. Interval type-2 fuzzy logic systems made simple. IEEE Transactions on Fuzzy Systems, 2006, 14(6):808-821 doi: 10.1109/TFUZZ.2006.879986 [30] Mo H, Wang F Y, Zhou M, Li R M, Xiao Z Q. Footprint of uncertainty for type-2 fuzzy sets. Information Sciences, 2014, 272:96-110 doi: 10.1016/j.ins.2014.02.092 期刊类型引用(9)
1. 吴琼，苏建徽，解宝，陈永辉，黄赵军，黄诚. 基于最大净功率输出的PEMFC阴极供气系统优化控制研究. 太阳能学报. 2024(02): 283-290 . 
百度学术2. 韩国鹏，詹天旸，戴朝华，刘卓，张杲，谢卓峰，陈维荣. 基于自适应滑模观测器的氢燃料电池空气供应系统泄露故障容错控制. 中国电机工程学报. 2024(05): 1848-1860 . 百度学术3. 陈家城. 燃料电池空气供给系统建模与控制. 机电技术. 2023(03): 70-73+91 . 百度学术4. 王伟，仲璐璐，刘洋，赵珺，王霖青. 基于模糊强化学习的空气压缩机群组协调预测控制. 控制与决策. 2023(08): 2183-2191 . 百度学术5. 吴庆，赵涛，佃松宜，郭锐，李胜川，方红帏，韩吉霞. 基于FPSO的电力巡检机器人的广义二型模糊逻辑控制. 自动化学报. 2022(06): 1482-1492 . 本站查看6. 赵涛岩，曹江涛，李平，冯琳，商瑀. 区间二型模糊免疫PID在环己烷无催化氧化温度控制系统中的应用. 化工学报. 2022(07): 3166-3173 . 百度学术7. 常新宇，富雪，任朝晖，王永富. PEMFC空气供给系统的智能建模与自适应控制. 控制工程. 2022(09): 1593-1599 . 百度学术8. 申奥勇，詹跃东. 基于二次型性能指标的燃料电池过氧比RBF-PID控制. 信息与控制. 2021(04): 498-504+512 . 百度学术9. 邱晓兰. 超混沌金融系统的模糊建模及其控制. 南阳理工学院学报. 2020(02): 37-41 . 百度学术其他类型引用(13)
-
下载:
计量
- 文章访问数: 2666
- HTML全文浏览量: 340
- PDF下载量: 822
- 被引次数: 22
