
在满足系统约束条件的前提下，需通过船舶能量管理系统对航行速度、发电机组启停状态和发电机组运行功率分配进行优化调度^{[6]}，其目标函数为
$${\rm{min}}C = {C_{\rm{e}}} + {C_{{\rm{pr}}}}$$ (1) 式中：C为船舶的总运行成本，其计量单位为货币单位monetary unit，简称m.u.；C_{e}为船舶电力系统的总运行成本；C_{pr}为船舶推进系统的总运行成本。

对于船舶电力系统，其总运行成本C_{e}为
$${C_{\rm{e}}}\! \!=\!\! \mathop \sum \limits_{j \!=\! 1}^T\! \mathop \sum \limits_{i\! =\! 1}^{{N_{\rm{e}}}} \left\{ {\Delta {T_j}\! \times \!{P_{ij}} \times \left[ {{t_{ij}} \!\times \!( {{F_{{\rm{C}}i}} \times {S_{{\rm{FC}}i}}( {{P_{ij}}} ) \!+\! {W_{{\rm{C}}ij}}} )} \right] \!+ \!{Q_{{\rm{C}}ij}}} \right\}$$ (2) 式中：j =1, 2, ···, T，为船舶航行时间，其中T为总航行时间，h；i =1, 2, ···, N _{e，}为发电机数量，其中N_{e}为发电机总台数；ΔT_{j}为第j个航行时间段，h；P_{ij}为第i台发电机ΔT_{j}时间段内产生的功率，MW；t_{ij}为计算系数，对于第i台发电机，若其运行则取1，否则取0；F_{Ci}为发电机消耗单位燃料的成本，m.u.；S_{FCi}为第i台发电机的特定燃料消耗曲线；W_{Cij}为ΔT_{j}时刻启动发电机的维护成本，m.u.；Q_{Cij}为ΔT_{j}时刻启动发电机的启动成本，m.u.。
其中，第i台发电机产生的功率P_{ij}为
$${P_{ij}} = {P_{{\rm{el, }}j}} \times {P_{n{\rm{, }}i}}/\sum\limits_{i' \in E} {P_{n{\rm{, }}i'}}$$ (3) 式中：P_{el, j}为已启动发电机在ΔT_{j}时间内的功率输出，MW；P_{n, i}为n时刻第i台发电机的输出功率，MW，其中n为机组运行时间；i'=1, 2, ···, E，为已启动的发电机数量，其中E 为已启动的发电机总台数；P_{n, i'}为第i'台已启动发电机的额定功率，MW。
发电机在单位时间内的燃料消耗F_{Ci}与输出功率P_{i}的关系方程为
$${F_{{\rm{C}}i}} = {a_{0i}} + {a_{1i}} \times {P_i} + {a_{2i}} \times P_i^2 + {a_{3i}} \times P_i^3$$ (4) 式中，a_{0i}，a_{1i}，a_{2i}，a_{3i}均为计算系数。
在船舶能量系统中，根据特定的燃料消耗S_{FCi}即可确定成本极值，以表示单位时间内输出单位功率所需的燃料：
$${S_{{\rm{FC}}i}} = {F_{{\rm{C}}i}}/{P_i}$$ (5) 对于推进系统，其总运行成本C_{pr}为
$$\begin{split} & {C_{{\rm{pr}}}} = \mathop \sum \limits_{j = 1}^T \mathop \sum \limits_{q = 1}^{{N_{{\rm{pr}}}}} \{ \Delta {T_j} \times {P_{qj}} \times [ {t_{qj}} \times ( {F_{{\rm{C}}q}} \times \\&\qquad {S_{{\rm{FC}}q}}( {{P_{qj}}} ) + {W_{{\rm{C}}qj}} ) ] + {Q_{{\rm{C}}qj}} \} \end{split}$$ (6) 式中：q=1, 2, ···, N_{pr}，为柴油机数量，其中N_{pr}为柴油机总台数；P_{qj}为第q台柴油机在ΔT_{j}时间内的输出功率，MW；t_{qj}为计算系数，对于第q台柴油机，若其运行则取1，否则取0；F_{Cq}为第q台柴油机在ΔT_{j}时刻内的燃料成本，m.u.；S_{FCq}为第q台柴油机的特定燃料消耗曲线；W_{Cqj}为在ΔT_{j}时刻启动柴油机的维护成本，m.u.；Q_{Cqj}为在ΔT_{j}时刻启动柴油机的启动成本，m.u.。
其中，第q台柴油机产生的功率P_{qj}为
$${P_{qj}} = {P_{{\rm{pr, }}j}} \times {P_{n{\rm{, }}q}}/\mathop \sum \limits_{q' \in Y} {P_{n{\rm{, }}q'}}$$ (7) 式中：P_{pr, j}为已启动柴油机在ΔT_{j}时间内的推力输出，MW；P_{n, q}为n时刻第q台柴油机的输出功率，MW；P_{n, q'}为第 q' 台已启动柴油机的额定功率，MW，其中q' =1, 2, ···, Y，Y为已启动的柴油机总台数。

1） 船舶发电机的功率约束。
$${P_{{\rm{min, }}i}} \leqslant {P_i} \leqslant {P_{{\rm{max, }}i}}{\rm{, }}\;\;\forall i$$ (8) 式中，P_{min, i}和P_{max, i}分别为第i台发电机的下限功率和上限功率。
2） 发电机组的突变功率约束。
$$\frac{{\left {{P_{ij}}  {P_{i{\rm{, }}j  1}}} \right}}{{\Delta {T_j}}} \leqslant {R_{{\rm{c}}i}}{\rm{, }}\;\;\;\forall i,j$$ (9) 式中：P_{i,j−1}为第i台发电机在第j−1个航行时间段内产生的功率；R_{ci}为第i台发电机的最大突变功率。
3） 电力负载和发电机组电力输出之间的平衡约束。
$$\mathop \sum \limits_{i = 1}^{{N_{\rm{e}}}} {t_{ij}} \times {P_{ij}} = {P_{{\rm{el, }}j}}{\rm{, }}\;\;\forall j$$ (10) 4） 发电机组的最小停机时间约束。
$$\left {{T_{{\rm{off}} \to {\rm{on, }}i}}  {T_{{\rm{on}} \to {\rm{off, }}i}}} \right \geqslant {T_{{\rm{off\_min, }}i}}{\rm{, }}\;\;\;\;\forall i$$ (11) 式中：
${T_{{\rm{off}} \to {\rm{on, }}i}}$ 为 第i台发电机的开机时刻；${T_{{\rm{on}} \to {\rm{off, }}i}}$ 为第i台发电机的关机时刻；${T_{{\rm{off\_min, }}i}}$ 为发电机的最小停机时间。5） 发电机组的最小连续运行时间约束。
$$\left {{T_{{\rm{on}} \to {\rm{off, }}i}}  {T_{{\rm{off}} \to {\rm{on, }}i}}} \right \geqslant {T_{{\rm{on\_min, }}i}}{\rm{, }}\;\;\;\;\;\forall i$$ (12) 式中，
${T_{{\rm{on\_min, }}i}}$ 为发电机的最小连续运行时间。6） 预防异常熄火约束。
当任意一台发电机突然停止工作时，其余发电机应能提供足够的电力输出。
$$\mathop \sum \limits_{i = 1}^{{N_{\rm{e}}}} {t_{ij}} \times {P_{{\rm{max, }}i}}  {\rm{max}}\left\{ {{t_{ij}} \times {P_{{\rm{max, }}i}}} \right\} \geqslant {P_{{\rm{el, }}j}}{\rm{, }}\begin{array}{*{20}{c}} {}&{\forall i,j} \end{array}$$ (13) 
1） 船舶柴油机的功率约束。
$${P_{{\rm{min, }}q}} \leqslant {P_q} \leqslant {P_{{\rm{max, }}q}}{\rm{, }}\;\;\;\forall q$$ (14) 式中：P_{min, q}和P_{max, q}分别为第q台柴油机的下限功率和上限功率；P_{q} 为第 q 台柴油机的输出功率。
2） 柴油机组的突变功率约束。
$$\frac{{\left {{P_{qj}}  {P_{q{\rm{, }}j  1}}} \right}}{{\Delta {T_j}}} \leqslant {R_{{\rm{c}}q}}{\rm{, }}\;\;\;\forall q,j$$ (15) 式中：P_{q, j−1}为第q台柴油机在第j−1个航行时间段内产生的功率；R_{cq}为第q台柴油机的最大突变功率。
3） 推进力负载与柴油机组动力输出之间的平衡约束。
$$\mathop \sum \limits_{q = 1}^{{N_{{\rm{pr}}}}} {t_{qj}} \times {P_{qj}} = {P_{{\rm{pr, }}j}}{\rm{, }}\;\;\;\;\forall j$$ (16) 4） 柴油机组的最小停机时间约束。
$$\left {{T_{{\rm{off}} \to {\rm{on, }}q}}  {T_{{\rm{on}} \to {\rm{off, }}q}}} \right \geqslant {T_{{\rm{off\_min, }}q}}{\rm{, }}\;\;\;\forall q$$ (17) 式中：
${T_{{\rm{off}} \to {\rm{on, }}q}}$ 为第q台柴油机的开机时刻；${T_{{\rm{on}} \to {\rm{off, }}q}}$ 为第q台柴油机的关机时刻；${T_{{\rm{off\_min, }}q}}$ 为柴油机的最小停机时间。5） 柴油机组的最小连续运行时间约束。
$$\left {{T_{{\rm{on}} \to {\rm{off, }}q}}  {T_{{\rm{off}} \to {\rm{on, }}q}}} \right \geqslant {T_{{\rm{on\_min, }}q}}{\rm{, }}\;\;\;\forall q$$ (18) 式中，
${T_{{\rm{on\_min, }}q}}$ 为柴油机的最小连续运行时间。 
船舶的有效推进力P_{epr}与航速V和船舶总阻力有关^{[7]}，其中总阻力的影响因素包括船舶自身属性和周围环境因素。通过改变船舶航速V，即可有效调整燃料消耗量，从而节约航行成本。船舶有效推进力P_{epr}与航速V的关系为
$${P_{{\rm{e}}  {\rm{pr}}}} = k{V^3}$$ (19) 式中，k为计算系数，本文取值为0.002 35。
在改变船舶航速V以减小燃油消耗的同时，航速应满足船舶的自身属性，即船舶航速的上限和下限约束条件：
$${V_{{\rm{min, }}j}} \leqslant {V_j} \leqslant {V_{{\rm{max, }}j}}{\rm{, }}\;\;\;\;\forall j$$ (20) 式中：V_{j}为船舶在ΔT_{j}时间段内的航速；V_{min, j}和V_{max, j}分别为船舶在ΔT_{j}时间段内的航速下限和上限。

1） 航行起点、终点的到达时间与距离约束。
假设船舶的出发时间为0时刻，到达终点港口的时间为T时刻，对于航行过程中的t时刻而言，其起始和终点距离约束为
$$t = 0,{L_0} = 0$$ (21) $$t = T = \mathop \sum \limits_{j = 1}^T \Delta {T_j},\;\;{L_T} = {d_{{\rm{total}}}}$$ (22) 式中：L_{0}为初始时刻的起点距离，n mile；j为初始时刻到终点时刻之间的时间离散点，即船舶航行时间；T为离散点总数，即总航行时间；ΔT_{j}为船舶航行过程中相邻2个离散点之间的时间间隔，即第j个航行时间段；L_{T}为终点时刻的航行总距离，n mile；d_{total}为实际航行总距离，n mile。
2） 航行至中间港口的到达时间与距离约束。
由于船舶并非从起点直接航行至终点，而是需要在航线所经港口装卸货物并补充食物和燃油，因此，船舶应在规定时间到达指定港口。
$${L_j} = {d_j}{\rm{, }} \;\;\;\;\forall j \in \tau $$ (23) 式中：L_{j}为到达各中间港口的总航行距离，n mile；d_{j}为船舶在各个中间港口之间的实际航行距离，n mile；τ为各个中间港口的时间离散点集合。

船舶能效运营指数（energy efficiency operational indicator，EEOI）即单位货物周转量所产生的温室气体排放量^{[8]}，代表了船舶在航行过程中的能耗效率，是评价其环境保护标准的一项重要指标。
在时间间隔为ΔT_{j}的船舶航行过程中，二氧化碳排放能效比
${\rho _{{\rm{EEO}}{{\rm{I}}_{1{\rm{, }}j}}}}$ 为$${\rho _{{\rm{EEO}}{{\rm{I}}_{1{\rm{, }}j}}}} = {m_{{\rm{C}}{{\rm{O}}_{2{\rm{, }}j}}}}/(\alpha \times {V_{{\rm{LF}}}} \times {V_j} \times \Delta {T_j})$$ (24) 式中：
${m_{{\rm{C}}{{\rm{O}}_{2{\rm{, }}j}}}}$ 为ΔT_{j}时间间隔内产生二氧化碳的总质量，g；α为船舶负载因数；V_{LF}为ΔT_{j}时间间隔内的货物运载量，t。船舶负载因数α与船舶类型密切相关，其计算公式为
$$\alpha = \left[ {\left( {0.1{n_{\rm{P}}}^\prime + {n_{\rm{V}}}^\prime } \right)/\left( {0.1{n_{\rm{p}}} + {n_{\rm{V}}}} \right)} \right]{G_{\rm{T}}}$$ (25) 式中：n_{P}'为船舶实际乘客人数；n_{V}'为船舶实际承载车辆数；n_{P}为船舶额定乘客人数；n_{V}为船舶额定承载车辆数；G_{T}为船舶总吨位，t。
由于二氧化碳排放总量与单位时间内柴油机和发电机的燃料消耗成本成正比^{[9]}，则有：
$$\begin{split} & {m_{{\rm{C}}{{\rm{O}}_{2{\rm{, }}j}}}} = \left\{ {\mathop \sum \limits_{i = 1}^{{N_{\rm{e}}}} {t_{ij}} \times {C_i} \times {S_{{\rm{FC}}i}}\left( {{P_{ij}}} \right) \times {P_{ij}} +}\right.\\&\;\left.{ \mathop \sum \limits_{q = 1}^{{N_{{\rm{pr}}}}} {t_{qj}} \times {C_q} \times {S_{{\rm{FC}}q}}\left( {{P_{qj}}} \right) \times {P_{qj}}} \right\} \times \Delta {T_j} \end{split}$$ (26) 式中：C_{i}为发电机的转换系数；C_{q}为柴油机的发电系数。
当船舶在港口停泊时，航速为0，则在时间间隔为ΔT_{j}的停泊过程中，二氧化碳排放能效比
${\rho _{{\rm{EEO}}{{\rm{I}}_{{\rm{2, }}j}}}}$ 为^{[10]}$${\rho _{{\rm{EEO}}{{\rm{I}}_{{\rm{2, }}j}}}} = {m_{{\rm{C}}{{\rm{O}}_{2{\rm{, }}j}}}}/(\alpha \times {V_{{\rm{LF}}}} \times \Delta {T_j})$$ (27) 根据船舶排放标准，船舶在航行及停泊时的二氧化碳排放能效比应小于环境要求，即
$${\rho _{{\rm{EEO}}{{\rm{I}}_{{\rm{1, }}j}}}} \leqslant {\rho _{{\rm{EEO}}{{\rm{I}}_{{\rm{1, max}}}}}}$$ (28) $${\rho _{{\rm{EEO}}{{\rm{I}}_{{\rm{2, }}j}}}} \leqslant {\rho _{{\rm{EEO}}{{\rm{I}}_{{\rm{2, max}}}}}}$$ (29) 式中，
${\rho _{{\rm{EEO}}{{\rm{I}}_{{\rm{1, max}}}}}}$ 和${\rho _{{\rm{EEO}}{{\rm{I}}_{{\rm{2, max}}}}}}$ 分别为船舶航行和停泊时二氧化碳排放能效比的环境保护标准上限值。 
船舶负载主要包括推进力负载和电力负载，且受船舶货物质量和乘客数量的影响。为了实现船舶柴油机和发电机的功率最优分配，需预测船舶的实时负载，从而准确控制柴油机和发电机的启停状态。对于给定的船舶电力负荷，可以制定发电机的启停及功率分配最优策略，并通过调整船舶航行速度来控制二氧化碳的排放量，从而满足航行及停泊时的环境保护要求。当实际负载与计算值不平衡时，采用船舶电力系统的实时最优控制方案^{[11]} 即可实现电力负载的再平衡，也可以采用负载优先启停策略^{[12]}对船舶能量系统进行调节。

粒子群优化算法（particle swarm optimizaton，PSO）是一种模仿鸟类在寻找食物过程中聚集情况的智能优化算法^{[13]}。鸟类在寻找食物的飞行途中会彼此传递信息，因此鸟类具备分散和聚集2种属性，虽然每只鸟的位置相对独立，但整体趋势却是向食物所在地不断聚集。在发现食物的过程中，虽然不知道食物的准确位置，但可以知道自身与食物的相对距离。通过比较每只鸟与食物之间的距离，即可得到离食物最近的一只鸟的位置，即该时间段的鸟群最优位置，也称之为粒子群的群体最优。通过信息分享，鸟群将以某种方式向群体最优的方向进行若干次靠拢，最终使绝大多数鸟类可以靠近食物所在位置，这实际上是所有可行解向最优解靠近寻优的过程。
优化问题目标函数f的最小值为
$$\min f\left( {{{{X}}_1}{\rm{, }}{{{X}}_2}{\rm{, }}{{{X}}_3}{\rm{, }} \cdots {\rm{, }}{{{X}}_J}{\rm{, }} \cdots {\rm{, }}{{{X}}_D}} \right){\rm{, }}\;\forall J \!=\! 1{\rm{, }}2{\rm{, }}3{\rm{, }} \cdots {\rm{, }}D$$ (30) 式中：X_{J}为优化向量的第J维分量；X_{D}为优化向量的第D维分量；J为1到D的整数，其中D为每个可行解的维数。
$${{X}}_J^L \leqslant {{{X}}_J} \leqslant {{X}}_J^U$$ (31) 式中，
${{X}}_J^L$ 和${{X}}_J^U$ 分别为每个可行解内某一维的下限和上限。 
1） 位置初始化。
粒子群优化算法是基于实数的编码寻优过程，首先需在可行解空间初始化，随机生成N_{P}个维度为D的可行解种群X_{I}，即
$${{{X}}_I} = \left( {{{{X}}_{I{\rm{,1}}}}{\rm{, }}{{{X}}_{I{\rm{,2}}}}{\rm{, }} \cdots {\rm{, }}{{{X}}_{I{\rm{,}}J}}{\rm{, }} \cdots {\rm{, }}{{{X}}_{I{\rm{,}}D}}} \right){\rm{, }}\; \forall I \!=\! 1{\rm{, }}2{\rm{, }}3{\rm{, }} \cdots {\rm{, }}{N_{\rm{P}}}$$ (32) 式中：I为迭代范围允许内的代数；X_{I, J}为第I代种群位置的第J维分量；X_{I, D}为第I 代种群位置的第D维分量；N_{P}为种群规模。由于初始种群在解空间内随机生成，故其可覆盖解空间内的所有范围。
2） 速度初始化。
完成可行解位置的初始化之后，即可对每个可行解的速度进行初始化：
$${{{V}}_I} = \left( {{{{V}}_{I,{\rm{1}}}}{\rm{, }}{{{V}}_{I,2}}{\rm{, }} \cdots {\rm{, }}{{{V}}_{I,J}}{\rm{, }} \cdots {\rm{, }}{{{V}}_{I,D}}} \right){\rm{, }} \;{\forall I \!=\! 1{\rm{, }}2{\rm{, }}3{\rm{, }} \cdots {\rm{, }}{N_{\rm{P}}}} $$ (33) 式中：V_{I}为第I 代种群的速度向量；V_{I, J}为第I 代种群速度的第J维分量；V_{I, D}为第I 代种群位置的第D维分量。
3） 个体最优和群体最优初始化。
完成可行解位置和速度的初始化之后，需要计算每个可行解的适应度值，即相当于鸟群中每只鸟与食物之间的距离，通过远近对比，即可得到距离食物最近的鸟的位置。因此，通过比较每个可行解的适应度值，即可得到种群初始化时的个体最优P_{best}和群体最优Z_{best}。
4） 速度位置更新。
得到可行解初始位置、初始速度、初始最优个体、初始最优群体值之后，即可利用粒子群公式对每个解的位置和速度进行刷新，则
$$\begin{split} & {{{V}}_{I + {\rm{1}}}} = w{{{V}}_I} + {c_1}ran{d_1}\left( {{{{Z}}_{{\rm{best}}I}}  {{{X}}_I}} \right) +\\&\qquad\;\;\;\; {c_2}ran{d_2}\left( {{{{P}}_{{\rm{best}}I}}  {{{X}}_I}} \right) \end{split}$$ (34) $${{{X}}_{I + {\rm{1}}}} = {{{X}}_I} + {{{V}}_I}$$ (35) 式中：V_{I+1}为第I +1代种群的速度向量；w为惯性系数；c_{1}和c_{2}为学习系数；rand_{1}和rand_{2}为介于[0，1]之间的随机数；Z_{bestI}和P_{bestI}
分别为第I代种群的群体最优值和个体最优值；X_{I+1}为第I +1代种群的位置向量。 可行解的速度和位置取值范围为：
$${{{V}}_{{\rm{min}}}} \leqslant {{{V}}_I} \leqslant {{{V}}_{{\rm{max}}}}$$ (36) $${{{X}}_{{\rm{min}}}} \leqslant {{{X}}_I} \leqslant {{{X}}_{\max}}$$ (37) 式中：V_{min}和V_{max}分别为粒子速度的下限和上限；X_{min}和X_{max}分别为粒子位置的下限和上限。
5） 适应度值更新。
根据每个粒子位置和速度的更新值，基于适应度值计算公式即可得到每个粒子适应度值的更新值。
6） 个体最优和群体最优更新。
将每一代更新之后的适应度值与个体最优和群体最优进行比较，如果当前适应度值小于此前个体最优和群体最优的适应度值，则以当前值取代此前的适应度值；否则，仍保留此前的适应度值。通过对比每一代的适应度值，即可刷新每一代个体最优和群体最优的位置，及其对应的每一代最优适应度值，从而实现个体最优和群体最优的更新。
7） 判断是否达到优化要求。
如果迭代次数达到最大迭代次数或优化结果满足输出条件，则迭代结束；否则，继续跳转至步骤4）进行迭代。

粒子群算法的权重系数将在很大程度上影响粒子群的迭代效率和搜索范围^{[14]}。标准粒子群算法一般采用固定的学习系数，或采用从大到小的线性迭代系数，但该算法不能较好地平衡全局最优和局部最优（容易陷入局部最优），且迭代速度较慢。因此，本文将采用非线性的学习系数，即初始权值较大，随着迭代次数的增加，其值将迅速减小。该自适应动态惯性权重w_{I}的计算公式为
$${w_I} = {{\rm{e}}^{  {\alpha ^I}/{\alpha ^{I  1}}}}$$ (38) 其中，
$${\alpha ^I} = \frac{1}{D}\sum\limits_{J = 1}^D {\left {f\left( {{{X}}_J^I} \right)  f({{X}}_{{\rm{min}}}^I)} \right} $$ (39) $${{X}}_J^I{\rm{ = }}{{X}}_{J{\rm{,1}}}^I,{{X}}_{J{\rm{,}}2}^I{\rm{, }} \cdots {\rm{, }}{{X}}_{J{\rm{,}}D}^I$$ (40) $$f({{X}}_{{\rm{min}}}^I) = \mathop {\min }\limits_{} f({{X}}_J^I)$$ (41) 式中：
${\alpha ^I}$ 和${\alpha ^{I  1}}$ 分别为第I代和I−1代种群适应度差值的均值；$f\left( {{{X}}_J^I} \right)$ 为第I代种群第J个解的适应度值，其中${{X}}_J^I$ 为第I代种群第J个解的向量，${{X}}_{J{\rm{,}}D}^I$ 为第I 代种群第J个解的第D维分量；$f({{X}}_{{\rm{min}}}^I)$ 为第I代种群所有可行解中的全局最优值，即第I代的Z_{best}。对w_{I}进行初始化时，令第1代w_{0}=0.729，则α^{−1}=3.16α^{0}。通过每一代所有可行解的适应度值与群体最优适应度值作差求平均值，即可判断函数的平整度。如果式（39）求得的α^{I}值较小，说明粒子的平整度较好；反之，则说明粒子的平整度较差。每一代都会产生新的适应度值，根据式（41）即可更新群体最优的适应度值，然后根据式（39）即可更新α^{I}值，从而得到自适应动态惯性权重参数。相较于线性递减的惯性权重，本文提出的惯性权重计算规则可以有效利用当前的收敛信息来选择合适的参数，这是一种带反馈的调节方式，增强了参数选择的指向性，且相对灵活准确。
鉴于不同迭代次数中α^{I}与α^{I−1}的比值可能波动较大，所以本文将其作为e的指数进行计算，从而保证每一代惯性权重参数w_{I}的平滑性，同时将其值限定于[0，1]范围内。当α^{I}/α^{I−1}>1时，说明这一代粒子整体趋于发散，比值越大，则发散程度越高，w_{I}越趋近于0，搜索步长越小；反之，当α^{I}/α^{I−1}<1时，说明这一代粒子整体趋于收敛，比值越小，则发散程度越低，w_{I}越趋近于1，搜索步长越大。通过比较判断每一代与上一代收敛发散程度的信息采集方式，即可明确每一代的惯性权重选择情况，从而实现更好的适应程度和优化效果。

步骤1：在可行解空间内初始化，随机生成N_{P}个D维的可行解种群，维度D包括柴油机和发电机的启停状态、船舶航速、发电机运行功率、柴油机发电功率等11维实数参数。
步骤2：生成种群初始速度和初始参数，包括最大迭代次数为G_{max}（其中初始代数G=1）、惯性系数w、学习系数c_{1}和c_{2}、随机数rand_{1}和rand_{2}等。
步骤3：对初始化参数进行调整与约束，包括发电机和柴油机的最大/最小功率约束、爬坡约束、最大/最小启停状态约束、EEOI约束等。
步骤4：计算初始适应度值，得出初始最优个体值和最优种群值。
步骤5：判断是否满足输出条件，即是否达到最大迭代次数或输出值满足期望值误差。如果满足，则输出结果；否则，进入步骤6。
步骤6：根据式（34）和式（35）进行速度、位置及适应度值更新。
步骤7：将每一代更新之后的适应度值与个体最优和群体最优进行比较，如果当前适应度值小于此前个体最优和群体最优的适应度值，则以当前最新的适应度值进行取代；否则，仍保留此前的适应度值。
步骤8：I=I+1，返回步骤5。

根据现有的研究成果，动态规划策略^{[5]}和差分进化（differential evolution，DE）算法^{[6]}都适用于船舶电力系统优化。本文将采用改进粒子群优化算法来调度船舶电力系统，并与文献[5]和文献[6]中的优化结果进行对比，用以验证本文算法的可行性与优化效果。文献[5]中：客渡船的总吨位为48 750 t，最高航速为23.5 kn，最大载客量为1 800人次，最大车载量为500辆，航行过程中二氧化碳排放量约束EEOI_{1}=21 g/(t·kn)，船舶进港停泊时二氧化碳排放量约束EEOI_{2}=120 g/(t·h)。为了与文献[5]中动态规划策略的优化结果进行对比，本文将采用相同的参数，包括船舶自身相关参数、乘客量及车辆数量、航行距离、船舶与港口之间的距离等。
在电力输出和推动力输出方面，本文为该船舶配置了3台发电机和2台柴油机，每台柴油机都通过轴连接至1台推动装置，具体参数如表1所示^{[5]}，其中P为发电机或柴油机的单位时间输出功率。
参数 发电机1 发电机2 发电机3 柴油机1 柴油机2 额定功率/MW 4 4 4 17.5 17.5 最大功率/MW 4 4 4 17.5 17.5 最小功率/MW 1 1 1 4.35 4.35 每克燃料释放的CO_{2}质量/g 2.5 2.5 2.5 3.2 3.2 最短启动时间和停止时间/h 1/1 1/1 1/1 1/1 1/1 系统启动和停止的消耗成本/m.u. 200/0 200/0 200/0 200/0 200/0 燃料消耗/（m.u.·t^{−1}） 500 500 500 450 450 每小时每单位功率的燃料质量
/kg·MW^{−1}·h^{−1}）343.5−80.3P+
12.5P^{2}346.7−73.8P+
11.12P^{2}345.6−69.6P+
10.44P^{2}211.7−5.21P+
0.231 5P^{2}210.7−4.47P+
0.198 8P^{2}Table 1. Data of ship power system
本文设定模拟航线的总长度为307.744 2 n mile，航线中间会停经3个港口，航行过程中船舶的实际载客数、实际货物量和负载因数如表2所示。
起止航线 实际载客数/人 实际载货量/t 负载因数 出发港—港口1 1 515 400 38 104 港口1—港口2 1 255 375 35 882 港口2—港口3 1 319 402 38 276 港口3—目的港 1 331 405 38 577 Table 2. Data for ship fullness
本文拟采用3种优化方案，具体如下：
1） 方案1：船舶配置3台发电机和2台柴油机，不采用最佳功率调度和推动力优化分配，不对EEOI进行限制约束。对于电力负载和推动力负载，发电机和柴油机都采用均匀分配的方式。
2） 方案2：船舶配置3台发电机和2台柴油机，采用最佳功率调度，不采用推动力优化分配，不对EEOI进行限制约束。
3） 方案3：船舶配置3台发电机和2台柴油机，采用最佳功率调度和推动力优化分配，对EEOI进行限制约束。

对于给定的船舶总负载、电力负载、推动力负载和航速，即可绘制船舶相关参数的曲线图，如图1所示。设定船舶的航行时间分为36个时间段（ΔT_{1}~
ΔT_{36}，下文简称T_{1}~T_{36}，图中则直接采用数字进行标示），每个时间段的间隔均为0.5 h。 方案1和方案2设定了航速不变，因此相应的船舶推动力负载也保持不变；在环保限制EEOI方面，方案1和方案2均未采取相关约束进行优化。但方案2的船舶发电机采用了电功率优化分配策略，可以在电力负载已知的情况下基于改进粒子群优化算法进行功率预分配，从而降低电力运行成本，而方案1仍采用了均匀分配电力负载的传统方式。方案3在方案2的基础上，进一步采用改进粒子群优化算法对船舶航速进行调整，用以优化船舶推进力负载，同时满足环保约束EEOI条件。
3种方案的优化对比情况如表3所示。
编号 电功率优化 推动力优化 航速调整 EEOI限制 1 × × × × 2 √ × × × 3 √ √ √ √ Table 3. Optimization comparison results of three schemes
图2是3种方案在T_{1}~T_{36}时间段的船舶航速对比示意图。由图2可知：未使用推动力优化策略的方案1和方案2的航速基本一致；方案3采用了推动力优化策略，其在到达港口和离开港口时的航速更高，而正常航行时的航速则明显降低，停泊时航速为0。
图3是3种方案在T_{1}~T_{36}时间段的船舶推进力负载对比示意图。由图3可知：船舶在T_{10}，T_{19}，T_{27}，T_{36}时刻的航速为0，推进力为0；与方案1和方案2相比，方案3在接近港口和离开港口时的推动力更高，而正常航行时的推动力负载则明显降低。
图4是3种方案在T_{1}~T_{36}时间段的船舶运行成本对比示意图。由图4可知：受船舶发电机和柴油机运行功率和启停工况的直接影响，燃料成本曲线基本与总负载正相关；与方案1和方案2相比，方案3经优化之后的总成本明显降低，这体现了改进粒子群优化算法的可行性。
图5是3种方案在T_{1}~T_{36}时间段的船舶环保约束EEOI值对比示意图。由于船舶在T_{10}，T_{19}，T_{27}，T_{36}时刻处于泊位，航速为0，所以采用EEOI_{2}限值，其余航行时刻则采用EEOI_{1}限值。为便于图形对比，将EEOI_{2}值缩小1/10进行处理。由图5可知，与方案1和方案2相比，方案3采取条件约束之后的限值可以满足环保标准要求。
图6所示为3种方案的3台发电机和2台柴油机在36个时间段内的启停状态及输出功率（假设所有发电机和柴油机已提前预热运行0.5 h）。
由图6（a）可知，对于未经优化的方案1而言，在电力负载（或推动力负载）低于单台发电机（或柴油机）最大功率的情况下，可以优先使用1台发电机（柴油机），例如由发电机1负责提供电力负载；反之，则均匀分配负载，即对每台启动的发电机（或柴油机）输出功率进行均分。在船舶从一个港口到达下一个港口的时间段内（T_{1}~T_{10}，T_{11}~T_{19}，T_{20}~T_{27}，T_{28}~T_{36}），船舶总负载呈现了中间高两端低的发展趋势，这与实船工况相符。
图6（b）中，方案2的推动力负载和总负载情况与方案1相似，但经电力负载优化分配之后，其整体表现更为平稳，可以在满足电力负载的前提下同时节约发电机的成本消耗。
由图6（c）可知，与方案1和方案2相比，方案3的推动力负载变化较大，通过推动力负载在各时间段的增减分配，最终使总负载基本维持在较稳定的水平。
根据仿真结果：方案1的运行成本为37 022.75 m.u.；采用电功率优化策略的方案2的运行成本为36 993.27 m.u.，降低了0.079 6%；采用电功率及推动力优化策略的方案3的运行成本为35 851.25 m.u.，降低了3.164 2%，且方案3可以满足EEOI环保要求。因此，改进粒子群优化算法可以降低船舶的运行成本并满足环保要求。

在相同的船舶参数和条件约束下，采用动态规划策略^{[5]}和差分进化算法^{[6]}时，方案2比方案1分别节约了0.05%和0.06%的成本，方案3比方案1分别节约了2.02%和2.90%成本；采用改进粒子群优化算法时，方案2比方案1节约了0.079 6%的成本，方案3比方案1节约了3.164 2%的成本。由此可见，在满足各项约束条件和环保要求的情况下，基于改进粒子群优化算法策略的船舶能量管理方案具有更好的经济性。
图7所示为改进粒子群优化算法与差分进化算法^{[6]}的运行成本收敛结果对比，可以看出：改进粒子群优化算法在迭代50次左右就基本实现了全局最优，且迭代效果更优，而差分进化算法则需要迭代150次，由此可见，改进粒子群优化算法具备更优的收敛性和稳定性。
Ship energy management scheme based on improved particle swarm optimization algorithm
doi: 10.19693/j.issn.16733185.01890
 Received Date: 20200226
 Rev Recd Date: 20200413
 Available Online: 20201110
 Publish Date: 20201230

Key words:
 marine power system /
 particle swarm optimization algorithm /
 greenhouse gas emission /
 energy management scheme
Abstract:
Citation:  YIN B, WANG X H, XIAO J M. Ship energy management scheme based on improved particle swarm optimization algorithm[J]. Chinese Journal of Ship Research, 2020, 15(6): 37–45 doi: 10.19693/j.issn.16733185.01890 