[发明专利]基于ALMBO优化算法的质子交换膜燃料电池子空间辨识方法

权利要求书

1.一种基于ALMBO优化算法的质子交换膜燃料电池子空间辨识方法，其特征在于，包括如下步骤：

步骤1，分析PEMFC工作特性和工作原理，选取合适的PEMFC状态空间辨识模型输入、输出变量；

步骤2，利用PEMFC测控平台采集PEMFC状态空间辨识模型建模所需的实验数据，并从采集的实验数据中挑选出有效实验数据，形成PEMFC状态空间辨识模型有效输入数据集与有效输出数据集；

步骤3，构造PEMFC状态空间辨识模型；

步骤4，PEMFC状态空间辨识模型最佳状态变量初始值寻优，具体描述为：

步骤4.1，初始化粒子：

生成NP组状态变量初始值x₀的随机值，记为NP是总的个体数，将中x₀的随机值平均分配到两个自定义种群Land1和Land2中，形成种群Land1的初始粒子x_o＝[s₁,s₂,…s_ceil(NP/2)]以及种群Land2的初始粒子x_v＝[s_ceil(NP/2)+1,s_ceil(NP/2)+2,…s_NP]，其中ceil为向上取整运算；

步骤4.2，变异反向学习迁移算子更新：

在经过变异反向学习时，将种群Land1中粒子x_o在第t+1时刻的更新值记为具体描述为：

式中，表示的第k维，即帝王蝶o的位置；lb_o、ub_o是搜索空间的下边界和上边界；p_r是变异概率，范围为(0.01，0.1)；R₁、b均为(0，1)间的随机数；表示在t时刻x_best的第k维，即帝王蝶在种群Land₁和Land₂当中最好个体的位置；t为当前迭代时刻；x_best为本次迭代得到的最优状态变量初始值；当b≤p_r时，通过随机反向学习扩大算法的搜索范围；当b＞p_r时，通过一般反向学习扩大算法的搜索范围；

步骤4.3，自适应调整算子更新：

完成粒子的迁移算子更新后，对Land2地中粒子x_v进行局部调整，将种群Land2中粒子x_v在第t+1时刻的更新值记为按式(3)～(5)进行更新：

当rand≤p时，

当rand＞p时，

其中，rand是(0，1)之间的一个随机数，表示的第k维，即帝王蝶v的位置；表示在t时x_best的第k维，即帝王蝶在种群Land₁和Land₂当中最好个体的位置；t为当前迭代时刻；x_best为本次迭代得到的最优状态变量初始值；表示从Land2中随机选取第t代帝王蝶x_r3的第k维；同时，在此条件下若rand＞BAR，则有：

其中：t为当前迭代时刻，为t时刻的粒子x_v；定义粒子x_v的随机步长为dx，dx_k则为粒子x_v第k维的随机步长，其步长为：α是权值，其值为：α＝S_max/t²，S_max为帝王蝶最大步长；BAR为调整率，其计算公式为：

BAR＝|mt-n| (6)

式中，t为当前迭代时刻，m、n是常数，其值计算方式为：

式中，t_m是最大迭代次数，BAR_max、BAR_min是调整率的上下界，其范围为[0,1]；

步骤4.4，柯西变异扰动更新：

将经步骤4.2～步骤4.3更新的粒子和按适应度值进行排序，将适应度排在最后的5个粒子按式(9)柯西变异扰动进行更新：

其中，cauchy(0,1)＝1/π(1+η²)，η为柯西随机变量生成函数：

η＝tan(π(ξ-0.5)) (10)

式中，ξ∈rand(0,1)，rand(0,1)是在区间[0,1]上服从均匀分布的随机数，表示在t时刻x_o的第k维；

更新完成后，将和中粒子按适应度值按从小到大进行排序，记为s^t+1；

步骤4.5，赋值计算：

将s^t+1中第一个粒子赋值给状态变量初始值x₀；

步骤5，基于有效输入数据集、有效输出数据集构造Hankel矩阵，具体描述为：

其中，Y_p为过去输出Hankel矩阵，U_p为过去输入Hankel矩阵，Y_f为未来输出Hankel矩阵，U_f为未来输入Hankel矩阵，p为过去Hankel矩阵的维数，f为未来Hankel矩阵的维数；

步骤6，求解PEMFC状态空间辨识模型的阶次，具体描述为：

步骤6.1，利用过去输出Hankel矩阵Y_p、过去输入Hankel矩阵U_p构造新矩阵W_p＝[Y_p U_p]^T；

步骤6.2，计算O_i＝(Y_f-Y_fU_f^T(U_fU_f^T)^-1U_f)(W_p-W_pU_f^T(U_fU_f^T)^-1U_f)W_p；

步骤6.3，令W₁＝I_h×h，对矩阵W₁O_iW₂进行SVD分解：W₁O_iW₂＝USV^T，其中，I表示为单位矩阵，h是Hankel矩阵的行数，j是Hankel矩阵的列数，U和V是酉矩阵，满足U^TU＝I，V^TV＝I；S是除去对角阵元素不为0，其余都为0的矩阵；

步骤6.4，根据S中对角线元素的大小来确定系统的阶次n：提取出S中对角线元素并从大到小排序，画出柱状图，选取数值处于中间的阶次作为该系统的阶次；

步骤7，将PEMFC状态空间辨识模型的辨识输出值与实际输出数据进行比较，若满足输出误差要求则直接输出，如不满足，则进行优化算法的继续迭代，具体描述为：

步骤7.1，计算广义能观矩阵Γ和Γ^⊥：

Γ＝W₁^-1U₁S₁^1/2 (11)

其中，U₁＝U(:,1:n)，指取步骤6.3中矩阵U的第1列到第n列，n是步骤6中确定的系统阶次；U₂＝U(:,n+1:end)，指取步骤6.3中矩阵U的第n列到最后一列所有数据；S₁＝S(1:n,1:n)，指取步骤6.3中S第1行到第n行中第1列到第n列的元素；

步骤7.2，计算状态空间辨识模型的系数矩阵A、C：

A＝ΓΓ^⊥ (12)

C＝Γ(1:l,:) (13)

式中，Γ^⊥表示矩阵Γ的正交矩阵，Γ(1:l,:)指取矩阵Γ的第1行到第l行所有数据，l为输出变量个数；

步骤7.3，构造M矩阵与L矩阵：

(L₁ L₂ … L_h)＝Γ^⊥ (15)

表示矩阵U_f的伪逆矩阵，是伪逆符号；

步骤7.4，计算状态空间辨识模型的系数矩阵B、D：

式中，h是Hankel矩阵的行数；

步骤7.5，将状态空间辨识模型的系数矩阵A、B、C、D以及状态变量初始值x₀代入状态空间模型，求解辨识输出值若达到预设最大迭代次数或者满足判决条件：

则输出辨识结果否则转步骤4，其中，J为输出误差。