多参数多目标混沌粒子群参数寻优方法

摘要

本发明公开了一种多参数多目标混沌粒子群参数寻优方法，步骤包括：步骤1，确定目标函数和待优化参数；步骤2，算法初始化；步骤3，计算种群中每个个体对应的目标函数；步骤4，更新个体历史最优解；步骤5，更新粒子速度和位置；步骤6，更新全局最优解集；步骤7，更新全局最优解；步骤8，结果判断。与一般的随机初始化方法和现有混沌Logistic映射粒子群初始化方法相比，本发明方法提高了全局寻优的性能且稳定性好；与一般的多目标加权的优化方法相比，本发明方法采用了Pareto最优解技术，解决了多目标加权方法中权重选取困难的问题。

说明书

技术领域

本发明属于优化控制技术领域，涉及一种多参数多目标混沌粒子群参数 寻优方法。

背景技术

近些年来，遗传算法、蚁群算法、免疫算法、粒子群优化算法等智能优 化算法在各个领域得到了广泛应用。其中粒子群优化算法，是通过模拟鸟群 觅食行为发展起来的一种基于群体协作的随机搜索算法，具有简单易行、收 敛速度快、优化效率高、鲁棒性好等特点，在处理优化问题中取得了很好的 效果，但是针对多参数、多目标综合寻优问题仍然是一个具有挑战性的课题。

目前，存在着两个主要问题，制约了现有粒子群方法在多参数、多目标 寻优问题中的应用。第一，针对多目标问题，现有的很多方法是将多目标乘 以权系数然后相加，将多目标优化问题变换成对多目标加权的单目标优化问 题，然而权值的选择困难使该方法难以达到很好的效果。第二，粒子群优化 算法中，初始粒子的分布对于优化算法的性能有着明显的影响，随着待优化 参数数量的增加，其目标解空间的维数也随之增加，非均匀分布的多维解空 间初始粒子，将会降低算法的全局收敛性能，使得算法收敛缓慢，并易陷入 局部最优解。

发明内容

本发明的目的是提供一种多参数多目标混沌粒子群参数寻优方法，针对 初始粒子在解空间分布影响粒子群算法收敛性的问题，通过采用单向耦合映 像格子时空混沌映射来初始化粒子的初始位置和速度，更好地实现粒子的初 始均匀分布；其次通过将粒子群算法与Pareto最优解理论相结合，得到多目 标函数最佳权衡基础上的解，解决了加权多目标方法中权值选取困难的问题。

本发明所采用的技术方案是，一种多参数多目标混沌粒子群参数寻优方 法，该方法按照以下步骤实施：

步骤1，确定目标函数和待优化参数

确定待优化参数数量为N，待优化目标函数个数ns，针对问题的不同优 化目标选取合适的目标函数J₁，J₂，…，J_ns；

步骤2，算法初始化

初始化粒子种群规模为M，粒子为多维粒子，其维数等于待优化参数数目 N，最大迭代次数为k_max，初始全局最优解集中解的个数j＝0，初始个体历史最 优解集中解的个数i_m＝0，随机初始化(第一次迭代时)全局粒子最优位置 Y_g＝[Y_g1,…,Y_gN]，初始化粒子的个体历史最优目标函数J_1max(m)＝0， J_2max(m)＝0，…，J_nsmax(m)＝0，m＝1,2,…,M，针对多个待优化参数，采用单向耦 合映像格子时空混沌模型初始化粒子群优化算法的粒子初始位置及初始速 度，初始位置表示为X_mn(0),初始速度表示V_mn(0)，其中m＝1,2,…,M，n＝1,2,…,N， 分别表示为第0次迭代时第m个粒子第n维的位置、飞行速度；

所述的单向耦合映像格子时空混沌映射模型为：

L_n(m+1)＝(1-ε_n)f[L_n(m)]+ε_nf[L_n-1(m)]，(1)

式中f[L_n(m)]为Logistic混沌映射，f[L_n(m)]＝μL_n(m)(1-L_n(m))，L_n(m)为状 态变量；n为空间格点位置，对应维数；m代表离散时间，对应种群规模；ε_n 为耦合强度，采用时空混沌映射，在优化算法开始阶段，直接生成具有多维 解空间均匀分布的初始粒子位置和速度，即选取粒子位置X_mn(0)＝k_x*L_n(m)， V_mn(0)＝k_v*L_n(m+M)，其中k_x和k_v分别为比例系数，将混沌时空模型的状态从[0， 1]间通过乘以对应的系数转换到对应位置和速度参数取值范围内；

步骤3，计算种群中每个个体对应的目标函数的值

将当前第k代粒子m,m＝1,…,M，所代表的多维参数X_mn(k)，n＝1,…,N，代 入待优化问题，计算该次迭代中粒子m对应的所有目标函数的值，表示为[J₁(m, k)，J₂(m,k)，…，J_ns(m,k)]；

步骤4，更新个体历史最优解

对于当前第k代中的第m个粒子，比较其对应的多个目标函数的值与其自 身历史最优位置对应的多个目标函数的值，如果该粒子的每个目标函数的值 都不劣于该粒子的历史最优位置对应的目标函数的值，那么用该粒子的目标 函数的值替换个体历史最好目标函数的值，同时将该粒子位置保存入个体最 优集中，具体过程是：

若J₁(m,k)≥J_1max(m)&J₂(m,k)≥J_2max(m)&…&J_ns(m,k)≥J_nsmax(m)，"&"表示 逻辑与操作，则J_1max(m)＝J₁(m,k)，J_2max(m)＝J₂(m,k)，…，J_nsmax(m)＝J_ns(m,k)， P_m(i_m，:)＝[X_m1(k),X_m2(k),…,X_mN(k)]，i_m＝i_m+1；其中P_m表示个体m的最优解集， i代表个体最优解集中个体最优解的数目，P_m(i_m，:)表示满足上述粒子m非劣 解条件从而保存入的个体最优解集中的第i_m个N维粒子；每个粒子都有一个个 体历史最优解集，其中保存了该粒子(在前k次迭代中)的历史最优解；

步骤5，更新粒子速度和位置

粒子速度和位置的更新迭代公式为：

V_mn(k+1)＝wV_mn(k)+c₁r₁(Y_m(n)-X_mn(k))+c₂r₂(Y_g(n)-X_mn(k))，(2)

X_mn(k+1)＝X_mn(k)+V_mn(k)，(3)

其中，X_mn(k)表示粒子m在第k代的第n维位置坐标；

V_mn(k)表示粒子m在第k代的第n维坐标上的速度(相当于位置变化率)；

c₁，c₂为学习因子，是非负常数；

r₁，r₂是介于[0，1]之间两个独立的随机数；w为惯性权重；

Y_m(n)从个体历史最优解集中随机选取，Y_m(n)＝P_m(ii_m，n)，ii_m＝random(1， i_m-1)，Y_g(n)为粒子第n维的全局最优解位置；

步骤6，更新全局最优解集

比较当前(第k代)种群中所有M个粒子对应的每个目标函数的值，将种 群中所有粒子中，任意一个目标函数为最大的粒子保存入Pareto全局最优解集 P_g中，同时该粒子对应的目标函数的值保存入Pareto最优目标函数值集合P_J中， 具体算法为：

若J₁(m,k)＝max{J₁(1,k),…,J₁(M,k)}orJ₂(m,k)＝max{J₂(1,k),…,J₂(M,k)} or…orJ_ns(m,k)＝max{J_ns(1,k),…,J_ns(M,k)}，其中“or”为逻辑或操作，则粒 子m的目标函数的值保存入Pareto最优目标函数值集中，即P_J(j，:)＝[J₁(m, k),…,J_ns(m,k)]，粒子m的位置保存入Pareto全局最优解集P_g(j，:)＝[X_m1(k),…, X_mN(k)]，j＝j+1，j代表Pareto最优解集中个体最优解的数目；

步骤7，更新全局最优解Y_g

求取Pareto全局最优解集中每个粒子对应的多个目标函数与其解集中各 目标函数最大值间的欧式距离，选取欧式距离最小的粒子作为全局最优解， 具体过程是：

选取全局最优解Y_g＝[Y_g1,Y_g2,…,Y_gN]＝P_g(jj,:)，其中jj为满足下式$ϵ\left(P_J\right(jj,:),$$[J_{g1},\mathrm{...},J_{gns}]=\underset{q=1,\mathrm{...},j-1}{\min }ϵ\left(P_J\right(q,:),[J_{g1},\mathrm{...},J_{gns}\left]\right)$的系数，其中[J_g1,…,J_gns]为 Pareto最优目标函数集中各目标函数分别能达到的最大值，即$J_{g1}=\underset{q=1,\mathrm{...}j-1}{\max }\{P_J$$(q,1)),J_{g2}=\underset{q=1,\mathrm{...}j-1}{\max }\{P_J(q,2)\},...,J_{{\mathrm{gns}}^{.}}=\underset{q=1,\mathrm{...}j-1}{\max }\{P_J(q,ns)\},ϵ(P_J(q,:),[J_{g1},...,J_{gns}])$表 示向量的欧式距离，P_J(q,:)表示Pareto最优目标函数集第q行所有元素；

步骤8，结果判断

k＝k+1，若达到停止条件，即k＝k_max，则返回最优粒子Y_g作为待优化问题 的解，同时返回其对应的目标函数的值；否则k＝k+1，返回步骤3重新更新。

本发明的有益效果是，采用单向耦合映像格子时空混沌模型初始化多维 粒子位置和速度，同时将粒子群优化算法与Pareto理论相结合，用以权衡多 目标寻优中有可能出现的多个目标解的优劣，以解决多目标之间的自动均衡 问题。与一般的随机初始化方法和现有混沌(Logistic映射)粒子群初始化方 法相比，本发明方法提高了全局寻优的性能且稳定性好；与一般的多目标加 权的优化方法相比，本发明方法采用了Pareto最优解技术，解决了多目标加 权方法中权重选取困难的问题。本发明方法不但能够应用在控制器参数寻优 中，同时可以应用在其他控制方法参数整定、电路拓扑参数优化等领域，具 有广泛的应用前景。

附图说明

图1是本发明方法采用的单向耦合映像格子时空混沌模型的示意图；

图2是本发明方法应用对象的复合有源箝位软开关三相功率因数校正变 换器示意图；

图3是本发明方法实施例中变换器双闭环PI控制器原理图；

图4是本发明方法实施例中变换器的目标函数计算结果示意图；

图5是本发明方法的流程图；

图6是现有方法针对多维粒子(N＝50时)不同粒子m所对应用于描述 分布均匀性的箱形图；

图7是本发明方法针对多维粒子(N＝50时)不同粒子m所对应用于描 述分布均匀性的箱形图；

图8为现有方法与本发明在选取粒子m＝30时N维元素的分布比较；

图9是采用不同优化方法所得PI参数控制变换器输出得到的直流电压波 形；

图10是采用不同优化方法所得PI参数控制变换器得到的输入A相电压、 电流波形；

图11是采用不同优化方法所得PI控制参数控制变换器在负载突变时的 输出直流电压波形；

图12是采用本发明方法优化所得PI参数控制变换器的输入A相电压、 电流波形实验结果；

图13是采用本发明方法优化所得PI参数控制变换器得到的负载突变时 的实验波形；

图14是采用现有PI参数整定方法优化所得PI参数控制变换器在负载突 变时得到的实验波形；

图15是采用加权混沌粒子群方法优化所得PI参数控制变换器在负载突 变时得到的实验波形。

具体实施方式

下面结合附图及具体实施例对本发明进行进一步的说明。

本发明针对初始粒子在解空间分布影响粒子群算法收敛性的问题，通过 采用单向耦合映像格子时空混沌映射来初始化粒子的初始位置和速度，更好 地实现粒子的初始均匀分布，提升了粒子搜索速度并避免陷入局部最优解。 同时将粒子群优化算法与Pareto理论相结合，用以权衡多目标寻优中有可能出 现的多个目标解的优劣，以实现多目标之间的自动均衡。

本发明方法作为一种针对多参数多目标情况下的混沌粒子群参数寻优的 方法，在此以复合有源箝位软开关三相功率因数校正变换器(以下简称变换 器)双闭环比例积分(PI控制器)为实施例对象，通过本发明参数优化方法 优化PI控制器参数以提高变换器的动静态性能。

下面对实施例对象(变换器及PI控制器)设计过程进行简要描述与分析：

参照图1，在复合有源箝位软开关三相功率因数校正变换器电路中，U_a、 U_b和U_c为三相输入电源电压，i_a、i_b和i_c为三相输入电流，L_a＝L_b＝L_c＝L为三 相输入交流滤波电感，三相输入电感及开关元件的等效电阻均为R，输出滤波 电容为C，直流侧负载为R_L；三相主电路桥臂采用6个IGBT元件(S₁-S₆) 两两串联后再并联而成，每个IGBT元件与一个二极管(D₁-D₆)对应两两反 并联，即每个IGBT元件的发射极与配对的二极管正极相联，每个IGBT元件 的集电极与配对的二极管负极相联；谐振电感L_r、辅助开关S₇和开关并联电 容C₁-C₇组成谐振电路，为软开关创造条件，C_c为箝位电容，箝位电路的加入， 降低了主开关两端的电压应力；七个IGBT元件的驱动采用12扇区空间矢量 调制算法进行控制，实现输入三相交流电压的整流及功率因数校正的要求，

根据电路基本规律(基尔霍夫电压、电流定律)，变换器在abc三相静止 坐标系下的模型为：

$\left(\begin{array}{c}L\frac{{\mathrm{di}}_a}{dt}=U_a-{\mathrm{Ri}}_a-\frac{2S_a-S_b-S_c}{3}(U_{dc}+U_{Cc})\\ L\frac{{\mathrm{di}}_b}{dt}=U_b-{\mathrm{Ri}}_b-\frac{2S_b-S_a-S_c}{3}(U_{dc}+U_{Cc})\\ L\frac{{\mathrm{di}}_c}{dt}=U_c-{\mathrm{Ri}}_c-\frac{2S_c-S_a-S_b}{3}(U_{dc}+U_{Cc})\\ L\frac{{\mathrm{di}}_{dc}}{dt}=S_a{i}_a+S_b{i}_b+S_c{i}_c-i_L\end{array}\right),---\left(4\right)$

其中，S_a,S_b,S_c分别表示三相桥臂的开关函数，当S_i＝1，i＝a,b,c时，代表 i相上桥臂IGBT导通，下桥臂IGBT关断；当S_i＝0，i＝a,b,c时，代表i相上 桥臂IGBT关断，下桥臂IGBT导通；U_dc为输出直流电压，U_Cc为箝位电容两 端电压，i_L＝U_dc/R_L为负载电流。

由三相abc坐标系变换到两相αβ坐标系的恒功率变换如下：

$\left(\begin{array}{c}{X}_{\alpha }\\ X_{\beta }\end{array}\right)=T_{abc/\alpha \beta }\left(\begin{array}{c}{X}_a\\ X_b\\ X_c\end{array}\right),$其中的变换矩阵$T_{abc/\alpha \beta }=\frac{2}{3}\left(\begin{array}{ccc}1& -\frac{1}{2}& -\frac{1}{2}\\ 0& \frac{\sqrt{3}}{2}& -\frac{\sqrt{3}}{2}\end{array}\right);$

再由两相αβ坐标系变换到dq旋转坐标系如下：

$\left(\begin{array}{c}{X}_d\\ X_q\end{array}\right)=T_{\alpha \beta /dq}\left(\begin{array}{c}{X}_{\alpha }\\ X_{\beta }\end{array}\right),$其中的变换矩阵w＝2πf为输入 正弦电压角速度；

通过上述的两次变换，将式(4)由abc坐标系下变换到dq坐标系下得 到式(5)：

$\{\begin{array}{c}L\frac{{\mathrm{di}}_d}{dt}=U_d-{\mathrm{Ri}}_d+{\mathrm{wLi}}_q-S_d(U_{dc}+U_{Cc})\\ L\frac{{\mathrm{di}}_q}{dt}=U_q-{\mathrm{Ri}}_q+{\mathrm{wLi}}_d-S_q(U_{dc}+U_{Cc})\\ C\frac{{\mathrm{dU}}_{dc}}{dt}=\frac{3}{2}(S_d{i}_d+S_q{i}_q)-i_L\end{array},---\left(5\right)$

由于C_c<<C，使得U_Cc<<U_dc，因此忽略U_Cc得到式(6)：

$\left(\begin{array}{c}L\frac{{\mathrm{di}}_d}{dt}=U_d-{\mathrm{Ri}}_d+{\mathrm{wLi}}_q-S_d{U}_{dc}\\ L\frac{{\mathrm{di}}_q}{dt}=U_q-{\mathrm{Ri}}_q-{\mathrm{wLi}}_d-S_q{U}_{dc}\\ C\frac{{\mathrm{dU}}_{dc}}{dt}=\frac{3}{2}(S_d{i}_d+S_q{i}_q)-i_L\end{array}\right),---\left(6\right)$

其中，U_d和U_q分别为dq坐标系下的d轴和q轴电压，i_d和i_q分别为dq 坐标系下的电流有功和无功分量，S_d和S_q分别为开关函数在dq坐标系下的开 关分量。

由式(6)得到的电流环方程如下：

$\{\begin{array}{c}L\frac{{\mathrm{di}}_d}{dt}=U_d-{\mathrm{Ri}}_d+{\mathrm{wLi}}_q-U_{rd}\\ L\frac{{\mathrm{di}}_q}{dt}=U_q-{\mathrm{Ri}}_q-{\mathrm{wLi}}_d-U_{rq}\end{array},---\left(7\right)$

其中，U_rd和U_rq分别为d轴和q轴的电流控制量，U_rd＝S_dU_dc，U_rq＝S_qU_dc。

电流环控制器采用前馈解耦PI控制策略如下：

$\{\begin{array}{c}{U}_{rd}=(k_{idp}{e}_{id}+k_{idi}{\int }_0^t{e}_{id}dt)+{\mathrm{wLi}}_q-{\mathrm{Ri}}_d+U_d\\ U_{rq}=(k_{iqp}{e}_{iq}+k_{iqi}{\int }_0^t{e}_{iq}dt)-{\mathrm{wLi}}_d-{\mathrm{Ri}}_q+U_d\end{array},---\left(8\right)$

其中e_id和e_iq分别表示对应的d轴和q轴电流误差，e_id＝i_d-i_dref，e_iq＝i_q-i_qref，i_dref和i_qref分别表示d轴和q轴电流参考值，k_idp、k_idi、k_iqp、k_iqi均为d轴和q轴电流PI 控制参数，

电流环控制器分为d轴PI控制器和q轴PI控制器两部分，d轴电流i_d为有功电 流，影响着变换器输入电流的幅值；d轴电流i_q为无功电流，影响着变换器的 功率因数，为了达到单位功率因数，需要i_qref＝0，i_dref由电压环控制器给出；可 见，电流环控制器的主要作用就是使输入电流为正弦且与输入电压同相位， 即保持变换器功率因数为1，因此，评价电流环PI控制器的优劣是看能否尽快 地使变换器达到单位功率因数，这需要d轴和q轴PI控制参数共同作用以达到良 好的控制效果。

由式(6)得到的电压环方程如下：

$C\frac{{\mathrm{dU}}_{dc}}{dt}=\frac{2}{3}(S_d{i}_d+S_q{i}_q)-i_L,---\left(9\right)$

电压外环采用如下控制器时，实现输出直流电压的跟踪控制，见下式：

$i_{dref}=k_{vp}(U_{dcref}-U_{dc})+k_{vi}{\int }_0^t(U_{dcref}-U_{dc})dt,---\left(10\right)$

其中U_dcref为电压参考值，k_vp、k_vi分别为电压环PI控制器比例与积分系数，

电压环控制器的主要作用是对输出直流电压进行调节，通常情况下，比 例系数k_vp起着加快电压响应速度，减少过渡过程时间的作用，但k_vp过大将引 起系统超调过大；积分系数k_vi起着减少输出电压稳态误差的作用，同样积分系 数影响超调和收敛速度。因此，考虑电压环性能是应从两个方面综合考虑其 性能，即输出直流电压的响应速度和超调，以及稳态误差。

综上所述，针对变换器输出直流电压和输入三相电流的不同控制目标， 式(10)的电压外环控制器和式(8)的电流环控制器共同组成的变换器实施 例的电压外环，电流内环的双闭环PI控制器，参照图2所示，其中(U_dcref-U_dc) 作为电压环控制器的输入，通过比例积分环节得到电压环控制器的输出i_dref； i_dref作为电流环d轴控制器的输入，与采样的三相输入电流经过三相变换后得到 的有功电流i_d计算得到d轴电流误差e_id，e_id通过电流环d轴比例积分环节，再加 上解耦项(U_d+wLi_q-Ri_d)得到d轴电流环的输出控制量U_rd；同理，I_qref＝0作为 电流环q轴控制器的输入，与采样的三相输入电流经过三相变换后得到的无功 电流i_q计算得到q轴电流误差e_iq，通过电流环q轴比例积分环节，再通过三相解 耦环节(U_q-wLi_d-Ri_q)得到q轴电流环的输出控制量U_rq；可见，双闭环控制器 其共有6个控制PI控制参数，这些控制参数共同作用以达到良好的控制效果。

参照图5，结合前述的实施例对象描述，本发明的方法，按照以下步骤具 体实施：

步骤1，确定实施例目标函数和待优化参数

待优化的目标函数为ns个，待优化参数的数量为N个，针对问题的不同 优化目标选取合适的目标函数，比如J₁，J₂，…，J_ns；

对于实施例对象，待优化参数为双闭环PI控制器参数，待优化参数PI 控制器参数数量N＝6，分别为k_vp、k_vi、k_idp、k_idi、k_iqp及k_iqi；待优化的目标函 数为ns＝3个，如图3所示，针对控制器需要满足的3个目标是：1)电压外 环作为控制器要求具有较快的响应速度和较小的超调，对应图3分隔线左部 分阴影区域尽可能小，该目标函数定义为J_v1；2)要求电压环的稳态误差为0， 对应图3分隔线右部分阴影区域趋近于0，即电压环的稳态误差趋近于0，目 标函数定义为J_v2；；3)电流内环应快速平稳的控制输入电流跟踪电压，实现 单位功率因数，即变换器功率因数与1间的差值尽可能小，目标函数定义为 J_i，则有：

$\begin{array}{cc}{J}_{v1}=\frac{1}{\frac{1}{L_1}\underset{i=1}{\overset{L_1}{\Sigma }}\left|e_i\right|}& \left|e_i\right|\ge ϵ\end{array},---\left(11\right)$

$\begin{array}{cc}{J}_{v2}=\frac{1}{\frac{1}{L_2}\underset{i=L_1}{\overset{L_2}{\Sigma }}\left|e_i\right|}& \left|e_i\right|\ge ϵ\end{array},---\left(12\right)$

$J_i=\frac{1}{\frac{1}{L_3}\underset{j=1}{\overset{L_3}{\Sigma }}\left|e_i\right|},---\left(13\right)$

其中，L₁代表输出直流电压在瞬态阶段的数据长度；L₂代表输出直流电压 在稳态阶段的数据长度；L₃代表功率因数值的数据长度；e_i和e_j分别表示输出直 流电压误差和功率因数误差，e_i＝U_dcref(iT)-U_dc(iT)，e_j＝PF(jT)-1；T为采样间隔， PF(jT)表示在第j个采样时间段的功率因数值；ε表示稳态误差的限，用于区分 瞬态和稳态阶段，误差大于ε时计入目标函数J_v1，否则计入目标函数J_v2。

步骤2，算法初始化

初始化粒子种群规模为M，粒子为多维粒子，其维数等于待优化参数的 数目N，最大迭代次数为k_max，初始全局最优解集中解的个数j＝0，初始个体 历史最优解集中解的个数i_m＝0，随机初始化(第一次迭代时)初始化全局粒 子最优位置Y_g＝[Y_g1,…,Y_gN]，初始化粒子的个体历史最优目标函数J_1max(m)＝0， J_2max(m)＝0，…，J_nsmax(m)＝0，m＝1,2,…,M，

针对多个待优化参数，采用单向耦合映像格子时空混沌模型初始化粒子 群优化算法的粒子初始位置及初始速度，初始位置表示为X_mn(0)，初始速度表 示V_mn(0)，其中m＝1,2,…,M，n＝1,2,…,N，分别表示为第0次迭代时第m个粒 子n维的位置、飞行速度，

对于上述的实施例对象，初始化粒子种群规模为M＝100，粒子为多维粒 子，其维数等于待优化参数的数目为N＝6，第k代的粒子m表示为X_m1(k)、 X_m2(k)、X_m3(k)、X_m4(k)、X_m5(k)和X_m6(k)，代表双闭环PI控制器六个控制参数 k_vp、k_vi、k_idp、k_idi、k_iqp和k_iqi。最大迭代次数为k_max＝10，随机初始化(第一次 迭代时)全局粒子最优位置Y_g＝[Y_g1,Y_g2,Y_g3Y_g4,Y_g5,Y_g6]∈[0,5]，初始化个体 历史最优目标函数J_v1max(m)＝0，J_v2max(m)＝0，…，J_nsmax(m)＝0，m＝1,2,…,M。 针对多个待优化控制参数，采用单向耦合映像格子时空混沌模型(式(1)所示) 初始化粒子群优化算法的粒子初始位置及速度，混沌时空模型初始位置取[0， 1]间的随机数，单向耦合映像格子时空混沌状态图，如图4所示，其生成的初 始粒子通过乘以系数k_x和k_v转换到对应参数范围[0，k_x]及[0，k_v]之间，对于粒 子位置初始化k_x＝5，对于粒子速度初始化k_v＝1，单向耦合映像格子时空混沌序 列耦合强度ε_n＝0.85，初始位置表示为[X_m1(0),X_m2(0),X_m3(0),X_m4(0),X_m5(0), X_m6(0)],初始速度表示[V_m1(0),V_m2(0),V_m3(0),V_m4(0),V_m5(0),V_m6(0)]，其中 m＝1,2,…,M。

步骤3，计算种群中每个个体对应的目标函数的值

将当前种群(第k代)中粒子m,m＝1,…,M，所代表的多维参数X_mn(k)， n＝1,…,N，代入待优化问题，计算该次迭代中粒子m对应的所有目标函数的 值，表示为[J₁(m,k)，J₂(m,k)，…，J_ns(m,k)]；

对于实施例，建立复合有源箝位软开关三相功率因数校正变换器双闭环 控制仿真模型。仿真模型参数设置如下：三相输入相电压U_in＝220V；三相输 入电感标称值L＝20mH；电感和开关等效电阻R＝1Ω；开关管并联电容 C₁＝C₂＝…＝C₇＝2nF；谐振电感L_r＝100μH；箝位电容C_C＝40μF；输出滤波电容 C＝1500μF；负载电阻R_L＝300Ω；输出期望电压U_dcref＝600V；开关频率 f＝10kHz。将当前粒子的六个参数[X_m1(k),X_m2(k),X_m3(k),X_m4(k),X_m5(k),X_m6(k)] 代入复合有源箝位软开关三相功率因数校正变换器仿真模型，运行仿真程序， 采集输出直流电压，电流波形及对应的功率因数值，从而计算得到粒子m在第 k代时对应的目标函数J_v1(m,k)，J_v2(m,k)，J_i(m,k)的值；

步骤4，更新个体历史的最优解

对于当前第k代中的第m个粒子，比较其对应的多个目标函数的值与其 自身历史最优位置对应的多个目标函数的值，如果该粒子的每个目标函数的 值都不劣于该粒子的历史最优位置对应的目标函数的值，那么用该粒子的目 标函数的值替换个体历史最好目标函数的值，同时将该粒子位置保存入个体 最优集中，具体过程是：

若J₁(m,k)≥J_1max(m)&J₂(m,k)≥J_2max(m)&…&J_ns(m,k)≥J_nsmax(m)，"&"表示 逻辑与操作，则J_1max(m)＝J₁(m,k)，J_2max(m)＝J₂(m,k)，…，J_nsmax(m)＝J_ns(m,k)， P_m(i_m，:)＝[X_m1(k),X_m2(k),…,X_mN(k)]，i_m＝i_m+1；其中P_m表示个体m的最优解 集，i代表个体最优解集中个体最优解的数目，P_m(i_m，:)表示满足上述粒子m 非劣解条件从而保存入的个体最优解集中的第i_m个N维粒子；每个粒子都有 一个个体历史最优解集，其中保存了该粒子(在前k次迭代中)的历史最优 解；

对于本实施例，具体过程是：若J_v1(m,k)≥J_v1max(m)&J_v2(m,k)≥J_v2max(m)& J_i(m,k)≥J_imax(m)，则J_v1max(m)＝J_v1(m,k)，J_v2max(m)＝J_v2(m,k)，J_imax(m)＝J_i(m,k)， P_m(i_m，:)＝[X_m1(k),X_m2(k),…,X_m6(k)]，i_m＝i_m+1；其中P_m表示个体m的最优解集， i_m代表个体m最优解集中个体最优解的数目，P_m(i_m，:)表示满足上述粒子m 非劣解条件从而保存入的个体m最优解集中的第i_m个N维粒子；

步骤5，更新粒子速度和位置

粒子速度和位置的更新迭代公式为

V_mn(k+1)＝wV_mn(k)+c₁r₁(Y_m(n)-X_mn(k))+c₂r₂(Y_g(n)-X_mn(k))，(14)

X_mn(k+1)＝X_mn(k)+V_mn(k)，(15)

式中X_mn(k)表示粒子m在第k代的第n维位置坐标；V_mn(k)表示粒子m在第k 代的第n维坐标上的速度(位置变化率)；c₁和c₂为学习因子，是非负常数； r₁和r₂是介于[0，1]之间两个独立的随机数；w为惯性权重；Y_m(n)从个体历史 最优解集中随机选取，Y_m(n)＝P_m(ii_m，n)，ii_m＝random(1，i_m-1)，Y_g(n)为粒子第 n维的全局最优解位置；

对于本实施例，学习因子c₁＝2.5，c₂＝1.5；惯性权重w＝0.9；速度限幅V_max＝1， r₁和r₂是介于[0,1]之间两个独立的随机数；w为惯性权重；Y_m(n)从个体历史最 优解集中随机选取，Y_m(n)＝P_m(ii_m，n)，ii_m＝random(1,i_m-1)，Y_g(n)为粒子的全局 最优解第n维位置；

步骤6，更新全局最优解集

比较当前第k代种群中所有M个粒子对应的每个目标函数的值，将种群 中所有粒子中，任意一个目标函数的值为最大的粒子保存入Pareto最优解集 P_g中，同时该粒子对应目标函数的值保存入Pareto最优目标函数值集合P_J中， 具体算法为：

若J₁(m,k)＝max{J₁(1,k),…,J₁(M,k)}orJ₂(m,k)＝max{J₂(1,k),…,J₂(M,k)} or…orJ_ns(m,k)＝max{J_ns(1,k),…,J_ns(M,k)}，其中“or”为逻辑或操作，则粒 子m的目标函数的值保存入Pareto最优目标函数值集中，即P_J(j，:)＝[J₁(m, k),…,J_ns(m,k)]，粒子m的位置保存入Pareto全局最优解集P_g(j，:)＝[X_m1(k),…, X_mN(k)]，j＝j+1，j代表Pareto最优解集中个体最优解的数目；

对于本实施例，更新全局最优解集的具体算法为：若J_v1(m,k)＝max{J_v1(1, k),…,J_v1(M,k)}orJ_v2(m,k)＝max{J_v2(1,k),…,J_v2(M,k)}orJ_i(m,k)＝max{J_i(1, k),…,J_i(M,k)}；则P_g(j，:)＝[X_m1(k),…,X_m6(k)]，所对应的Pareto最优目标函 数集P_J(j，:)＝[J_v1(m,k),J_v2(m,k),J_i(m,k)]，j＝j+1，j代表Pareto最优解集中 个体最优解的数目；

步骤7，更新全局最优解Y_g

求取Pareto最优解集中每个粒子对应的多个目标函数与其解集中各目标 函数最大值间的欧式距离，选取欧式距离最小的粒子最为全局最优解，具体 过程是：

选取全局最优解Y_g＝[Y_g(1),Y_g(2),…,Y_g(N)]＝P_g(jj,:)，其中jj为满足下式$ϵ$$\left(P_J\right(jj,:),[J_{g1},...,J_{gns}\left]\right)=\underset{q=1,\mathrm{...}j-1}{\min }ϵ\left(P_J\right(q,:),[J_{g1},...,J_{gns}\left]\right)$的系数，其中[J_g1,…, J_gns]为Pareto最优目标函数集中各目标函数分别能达到的最大值，即 $J_{g1}=\underset{q=1,\mathrm{...}j-1}{\max }\left\{P_J(q,1)\right\},J_{g2}=\underset{q=1,\mathrm{...}j-1}{\max }\left\{P_J(q,2)\right\},...,J_{gns}=\underset{q=1,\mathrm{...}j-1}{\max }\left\{P_J(q,ns)\right\},ϵ\left(P_J\right(q,:),$$[J_{g1},...,J_{gns}])$表示向量的欧式距离，P_J(q,:)表示Pareto最优目标函数集第q行 所有元素；

对于本实施例，全局最优解Y_g的更新过程表达为：

选取全局最优解Y_g＝[Y_g(1),Y_g(2),…,Y_g(6)]＝P_yg(jj,:)，其中jj满足下式$ϵ$$\left(P_J\right(jj,:),[J_{gv1},J_{gv2},J_{gi}\left]\right)=\underset{q=1,\mathrm{...},j-1}{\min }ϵ\left(P_J\right(q,:),[J_{gv1},J_{gv2},J_{gi}\left]\right),$其中[J_gv1,J_gv2,J_gi] 为Pareto最优目标函数集中各目标函数分别能达到的最大值，即 $J_{gv1}=\underset{q=1,\mathrm{...}j-1}{\max }\left\{P_J(q,1)\right\},J_{gv2}=\underset{q=1,\mathrm{...}j-1}{\max }\left\{P_J(q,2)\right\},J_{gi}=\underset{q=1,\mathrm{...}j-1}{\max }\left\{P_J(q,3)\right\},ϵ\left(P_J\right(q,:),$$[J_{gv1},J_{gv2},J_{gi}])$表示向量的欧式距离，P_J(q，:)表示Pareto最优目标函数集第q 行所有元素；

步骤8，更新k＝k+1，若达到停止条件，即k＝k_max，则返回最优粒子Y_g作 为待优化问题的解，同时返回其对应的目标函数的值；否则k＝k+1，返回步骤 3，重新进行更新。

本发明多参数多目标混沌粒子群优化方法，能够获得满足多目标控制性能 综合要求的一组PI控制参数，该方法初始粒子能够在多维解空间中实现混沌 分布，其效率和稳定性要优于一般随机初始化以及现有混沌映射初始化方法， 同时解决了现有加权多目标优化中权值难以选择的问题，能够大大缩短参数 调制过程，降低对于设计者的经验要求。

验证本发明方法的效果

图6为现有Logistic混沌映射和本发明单向耦合映像格子时空混沌映射分 布比较图。图6和图7是两种方法针对多维粒子(N＝50时)不同粒子m所对 应用于描述分布均匀性的箱形图，图中箱格的长度越长表明该粒子的N维元 素分布越均匀。从图6与图7的比较中可以看到，对于每一个粒子m，本发 明所采用时空混沌模型粒子多维元素分布的均匀性都要优于现有Logistic混 沌映射。图8为两种方法当选取粒子m＝30时N维元素的分布情况，可以看 到，相对于现有Logistic映射，单向耦合映像格子时空混沌映射具有更为严格 的边界条件，即图8中采用Logistic方法产生的多维粒子包含的“星点”代表 的控制参数为0，实际控制中是不合理的，其在单向耦合映像格子时空混沌映 射中是没有出现的。

本发明所提方法具有产生粒子在多维解空间中均匀分布的优点，有助于 提升优化算法稳定性及加快算法收敛速度，为了证明这一点，基于同样的多 目标粒子群环境，采用本发明时空混沌模型粒子初始化方法，随机数初始化 方法和现有Logistic映射初始化方法的仿真比较结果如表1所示。表1结果为 将上述方法均多次运行(50次)后的结果比较，其中k_r表示不同方法达到优 化稳定解所需的平均迭代次数，可以看到，本发明方法具有最快的优化搜索 速度和最大平均目标函数值，即优化效果最好。

表1、不同粒子初始化方法效果比较

为了验证本发明方法基于Pareto最优解的多目标混沌粒子群的有效性， 与现有PI控制参数整定方法，加权多目标混沌粒子群方法(目标函数为J＝2J_v1+ 3J_v2+J_i)的仿真比较结果如表2所示。

表2、三种不同方法优化结果比较

从表2可以看到本发明方法优化所得控制参数目标函数的值最大，优化 效果最好。

采用表2不同优化方法所得PI控制参数的复合有源箝位软开关三相功率 因数校正变换器的输出直流电压波形仿真结果如图9所示，所对应的输入A 相电压、电流波形如图10所示。由图9和图10可见，采用本发明优化算法 优化所得的控制参数具有最好的控制效果，其不需依靠设计者经验不断调整 参数，也没有加权参数难以选择的问题，其输出直流电压具有最小的超调并 且没有稳态误差，输入A相电流更快的达到了单位功率因数。

负载由300欧姆突变到450欧姆时，采用不同优化方法所得PI控制参数 的仿真结果如图11所示，可以看到，采用本发明优化方法所得PI控制参数具 有最好的控制效果。

搭建复合有源箝位软开关三相功率因数校正变换器的样机，样机电路参 数与仿真一致，控制算法采用DSP28335数字控制器实现。得到实验结果如下： 采用本发明优化方法优化所得PI控制参数控制复合有源箝位软开关三相PFC 变换器得到的输入A相电压、电流波形如图12所示。图12为输入A相电压、 电流波形，其中横坐标为时间，单位为秒，A相电压的纵坐标单位为伏特，A 相电流的纵坐标单位为安培，可见输入电流波形为正弦，能够跟随输入电压 变化，即达到了单位功率因数，采用HIOKI3197功率因数分析仪测量得到， 三相平均功率因数为0.991。

本发明参数优化方法所得PI控制参数负载突变(R_L由300欧姆变至450 欧姆)的实验结果如图13所示，其中通道1(channel1)为输出直流电压波 形，横坐标为时间，单位为秒，纵坐标为电压，单位为伏；通道2(channel2) 为输入A相电流波形，横坐标为时间，单位为秒，纵坐标为电流，单位为安；

现有PI参数整定方法所得PI控制参数负载突变(R_L由300欧姆变至450 欧姆)的实验结果如图14所示，其中通道1(channel1)为输出直流电压波 形，横坐标为时间，单位为秒，纵坐标为电压，单位为伏；通道2(channel2) 为输入A相电流波形，横坐标为时间，单位为秒，纵坐标为电流，单位为安；

现有加权混沌粒子群方法所得PI控制参数负载突变(R_L由300欧姆变至 450欧姆)的实验结果如图15所示，其中通道1(channel1)为输出直流电压 波形，横坐标为时间，单位为秒，纵坐标为电压，单位为伏；通道2(channel 2)为输入A相电流波形，横坐标为时间，单位为秒，纵坐标为电流，单位为 安；

对比图13、图14和图15可见，本发明方法优化所得PI控制参数在负载 突变时，电压扰动更小，恢复更快，性能更好。