一种基于误差反馈的阻力系数标定方法与流程

本发明属于燃料电池领域，具体涉及一种基于误差反馈的阻力系数标定方法。

背景技术：

1、燃料电池系统由空气回路、氢气回路和冷却回路三个子系统构成。三个回路的电堆、中冷器、加湿器、气水分离器等部件（bop）共同运作使系统正常运行。但是，bop内部结构的几何形状极为复杂，且单个电堆往往包含了数百个重复单元，如果在计算流体动力学（cfd）仿真中对每一个几何形状都进行建模，计算量将不可接受。因此工程上的做法是采用多孔介质等效模型，引入粘性阻力系数和惯性阻力系数来描述气体通过该区域受到的阻力。描述流体在多孔介质中流动压降的标准方程为darcy-forchheimer方程：

2、 (1)

3、其中，为单位长度压降，；为惯性阻力系数，；为渗透率，，其倒数为粘性阻力系数，；为流体密度，；为流体动力粘度，；为流体流速，。

4、和为经验参数，在现有技术方案当中，需要由工程师在实验台上测量不同流量下气体通过被测多孔区域结构的压降数据，再代入上述方程进行反算。通过在计算流体力学软件中赋予和的初始值进行数值模拟，比较模拟值与实验值压降的偏差手动修改参数后重新提交仿真，循环往复，直到残差收敛到可接受范围。

5、现有方案通常采用人工手动迭代标定法，手动修改阻力系数的方式效率低下，单次数值模拟得出结果通常需要数十分钟乃至数小时，人工试算的方式需要迭代数十次，同时对和进行经验性调整往往伴随着修正方向的随机性，因此难以获得全局最优解。另外，人工手动迭代标定的方式以单一目标电流点作为验证基准，未说明在全工况范围内标定结果是否保持一致，具体参见cn 115911467 a。另一种方式是采用理想化的几何建模，采用均匀n×n圆柱阵列来表征多孔介质几何结构，通过对比不同积分边界下的总体阻力系数误差，筛选出不受积分边界影响的“稀疏型”多孔介质模态，从而确定合理的固体率范围及相邻圆柱间的间距，并据此完成模型构建，此模型完全由数值仿真推算，未纳入实验测量数据进行校正，会使得真实流阻与理论计算存在系统性偏差，具体参见cn 116956771 a。

技术实现思路

1、本发明克服上述缺陷，提供一种基于误差反馈的阻力系数标定方法。本发明结合历史差商割线法预测阻力系数修正步长，并采用算法全自动化标定计算流体动力学当中反复标定阻力系数的过程，有效抑制复杂流场计算时候的数值崩溃。

2、本发明技术方案如下。

3、一种基于误差反馈的阻力系数标定方法，包括以下步骤：

4、步骤一、对氢燃料电池系统进行台架测试，获取系统测试温度、操作压力、不同工况下阻力件的实验压降值和流量值。所述阻力件包括但不限于中冷器、增湿器、电堆或气水分离器等；所述不同工况下具体为：包括但不限于工质相对湿度、操作压力或温度。

5、进一步地，由于各bop部件结构形状不尽相同，布置的压力测点接近于符合系统测试的满足得出可靠的测量值的必要位置即可，满足测量流体主流方向的总压位置。

6、步骤二、在计算流体动力学数值求解系统（例如基于有限体积法的开源流体力学软件openfoam）中设置粘性阻力系数和惯性阻力系数，设置当前数值求解迭代步的边界条件，并对求解对象进行数值求解；动态计算当前工况下流体物性参数，提取计算收敛后的阻力件区域的仿真压降值、局部特征流速、粘性阻力系数和惯性阻力系数值；

7、步骤三、计算数值求解结果的阻力件仿真压降值与测试阻力件压降目标值的相对误差，判断相对误差是否满足预设的收敛条件；

8、所述相对误差计算方式为：

9、 (2)

10、判断相对误差是否；其中，为相对误差的给定阈值，为多孔区域实验压降目标值；为多孔区域压降仿真值；

11、如果满足则输出阻力件的仿真压降值、粘性阻力系数和惯性阻力系数，结束当前工况的寻优；如果相对误差不满足收敛条件，则判断当前迭代步数是否达到预设的最大迭代限制；

12、如满足迭代步数，则触发发散保护机制结束标定，输出报错信息并结束寻优；如满足迭代步数，则进入步骤四；

13、步骤四、引入无量纲局部阻力比提取所述阻力件的局部流速进行物理流态判定，分别执行独立的阻力系数修正过程：

14、判断当前流场处于粘性主导的低速流态时，进入低速修正过程，锁定惯性阻力系数，计算并更新下一代的粘性阻力系数修正步长；

15、判断当前流场处于惯性主导的高速流态时，进入高速修正过程，锁定粘性阻力系数，计算并更新下一代的惯性阻力系数修正步长；

16、根据所述的流体在多孔介质中流动压降的标准方程，darcy-forchheimer方程，压降由惯性项压降和粘性项压降组成，公式变换为：

17、 (3)

18、令，，则；

19、其中，为总阻力降，；为惯性阻力，由系数控制，在高流速下固体通过多孔介质骨架时容易发生撞击和产生漩涡，此时消耗的能量与速度的平方成正比，此时由这种能量耗散主导流动；为粘性阻力，由系数控制，在低流速下流体在多孔介质中运动时容易与壁面发生摩擦，此时由这种摩擦主导流动；为了避免流速影响从而模糊了粘性阻力和惯性阻力主导的流态，定义阻力比γ：

20、 (4)

21、其中，为无量纲阻力比。通过对无量纲阻力比的大小判断，算法开启阻力系数的修正过程：

22、当时，流体的压降将由粘性阻力系数主导，判定当前流态低速流态，执行低速修正过程，锁定，修改。

23、当当前迭代步数时，采用当前步的相对误差、当前步的阻力系数与预设的初始比例常数连乘计算当前步数的粘性阻力系数修正步长：

24、 (5)

25、当当前迭代步数时，采用基于历史迭代的差商割线算法计算当前步数的粘性阻力系数修正步长：

26、 (6)

27、当时，流体的压降将由惯性阻力系数主导，判定当前流态为高速流态，执行高速修正过程，锁定，修改；

28、当当前迭代步数时，采用当前步的相对误差、当前步的阻力系数与预设的初始比例常数连乘计算当前步数的惯性阻力系数修正步长：

29、 (7)

30、当当前迭代步数时，采用基于历史迭代的差商割线算法计算当前步数的惯性阻力系数修正步长：

31、 (8)

32、预设的初始比例常数k，取值范围为 (0, 1]；

33、步骤五、结合自适应松弛因子对步骤四的阻力系数修正步长平滑处理后，计算并更新新一代阻力系数，将其输入至数值求解系统中的多孔介质区域设置中，并返回步骤二，直到形成闭环寻优；

34、现在开始执行自适应松弛因子的计算：

35、当当前迭代步数时，跳过乘积判断，当前步松弛因子保持初始预设值不变；

36、当当前迭代步数时，计算当前步数的相对误差与上一迭代步数的相对误差的乘积：，同时预设两个静态超参数：作为加速因子的第一调节系数，该系数优选为1.2；以及作为减速惩罚的第二调节系数，该系数优选为0.5；

37、如果乘积，判断系统处于“单调逼近区间”，按照第一调节系数放大所述松弛因子，令当前松弛因子为：

38、 (9)

39、其中,为当前迭代步松弛因子，为上一迭代步松弛因子；

40、如果乘积，当前已越过寻优目标值，判断系统处于“过冲区间”，按照第二调节系数缩小所述松弛因子，另当前松弛因子为：

41、 (10)

42、将计算得到的强制截断在(0.3, 0.7)的安全区间内。

43、确定当前步数的修正松弛因子后，根据步骤五的选择结果执行新一代阻力系数的更新：

44、若当前处于低速修正过程，那么下一代阻力系数计算公式为：

45、 (11)

46、 (12)

47、若当前处于高速修正过程，那么下一代阻力系数计算公式为：

48、 (13)

49、 (14)

50、系数更新完成后，令迭代步数，将新一代阻力系数输入数值求解系统中，返回步骤二执行闭环迭代寻优。

51、上述方法中，步骤二中，所述初始化粘性阻力系数和惯性阻力系数包括人为经验确定或基于步骤一测试数据进行标定。

52、上述方法中，步骤三中，所述相对误差的给定阈值范围在[0,1]，工程上优选为0.05，即5%。

53、上述方法中，步骤四中，预设的初始比例常数k为0.1。

54、上述方法中，步骤五中，所述自适应松弛因子初始预设值为0.5。

55、与现有技术相比，本发明的优势在于：

56、（1）本发明结合历史差商割线法预测阻力系数修正步长，并采用算法全自动化标定计算流体动力学当中反复标定阻力系数的过程，有效抑制复杂流场计算时候的数值崩溃，将长周期的人工标定替代为高效率的算法自动化标定。

57、（2）本发明中引入无量纲阻力比进行流态判定，解决了传统标定过程中将无关系数放大，使得求解过程残差居高不下的弊端，实现了高稳定性的参数解耦寻优。