CN115017774A - 利用改进abc算法确定热学参数的大坝温度预测方法 - Google Patents

Abstract

本发明涉及利用改进ABC算法确定热学参数的大坝温度预测方法，包括：采集大坝不同部位的混凝土温度监测数据；建立各个混凝土浇筑仓的有限元模型；确定热学参数，并建立目标函数；利用有限元模型进行温度仿真模拟计算，并计算目标函数值；利用改进的人工蜂群算法对热学参数进行迭代反演，判断得到的热学参数值的优劣；重复反演计算，得到满足要求的热学参数的取值，并将计算得到的温度值与混凝土实际温度值相比较，达到预期目标后停止迭代反演；根据得到的热学参数值，预测大坝混凝土的温度。本发明利用改进ABC算法对热学参数进行反演计算，减小了计算的热学参数值与真实值的误差，提高了大坝温度预测的精度，有利于对大坝实施精细化的实时温度调控。

Description

利用改进ABC算法确定热学参数的大坝温度预测方法

技术领域

本发明属于水利水电混凝土大坝温度调控领域，具体涉及一种利用改进ABC算法确定热学参数的大坝温度预测方法。

背景技术

在混凝土大坝的建设过程中，由于自身水化热和外界气温等因素的影响，坝体中容易产生较大的温度应力，当混凝土的抗拉强度不足以抵抗该温度应力时，极易产生温度裂缝。为了实时准确地了解大坝的温度性态，仅靠现有仪器得到当前的温度数据是不够的，还需要通过软件仿真对大坝温度未来变化趋势进行预测，获得准确的混凝土温度场来分析大坝应力状态。然而，混凝土温度场仿真计算过程中易受多种不确定因素影响，如外部温度、内部水管冷却以及混凝土的热学参数等。外部温度、冷却水管流量可现场监测得到，而热学参数取值则需要通过室内试验获得。但由于室内试验的局限性，热学参数的室内实验值与实际值存在较大误差。因此，有必要基于大坝施工现场实测数据和智能优化分析方法，综合考虑现场外界气温与冷却通水因素的影响，通过反演计算获得大坝施工过程中真实热学参数值。

大坝热学参数反演计算过程实质上是利用优化算法调用温度场仿真计算程序进行多次迭代的过程。现有技术的优化算法调用温度场的方法局部搜索能力差、对初始值依懒性强，在复杂问题上，容易陷入局部最优、无法保证能够搜索到全局最优解，并且还存在计算量大、迭代时间长等问题。

人工蜂群算法(Artificial Bee Colony algorithm，ABC)是一种模拟蜜蜂群体寻找优质蜜源的仿生智能计算方法，由土耳其学者Karaboga在2005年提出。它是通过对自然界蜂群的长期观察，将蜂群的群体智能用于求解优化问题并与人工智技术相结合而产生的。该算法属于集群思想的范畴，类似的还有鱼群算法、蚁群算法等。人工蜂群算法本质是一种广义的邻域搜索算法，它将蜂群分为引领蜂、跟随蜂、侦察蜂3种蜜蜂类型以及分享蜜源和放弃蜜源2种基本行为。每个蜜源的位置代表优化问题的可行解，蜜源的花蜜量代表对应解的质量，蜜蜂寻找蜜源的速度代表求解对应的优化问题的速度，通过引领蜂、跟随蜂、侦察蜂在不同情况下的转换，借助启发式策略，不仅能有效进行局部搜索，还具有全局寻优的能力。

与其它智能算法相比，ABC算法具有控制参数少，计算简单易于实现，搜索速度较快，以适应度函数作为主要进化依据等优良特点，但是ABC依然存在收敛速度后期放慢以及易陷入局部最优解等问题。

因此，研究改进的人工蜂群算法，并将其应用到对大坝热学参数的反演计算以实现对大坝温度的精准预测。

发明内容

本发明的目的是针对上述问题，提供一种利用改进ABC算法确定热学参数的大坝温度预测方法，在ABC算法解空间的搜寻过程中增加当前全局最优值的引导作用以加快ABC算法的收敛速度，并在引领蜂领域搜索过程中增加交叉操作以提高全局寻优能力；利用改进ABC算法对热学参数进行反演计算，减小计算的热学参数值与真实值的误差；再根据计算得到的热学参数值预测大坝温度，以便于对大坝实施精细化温度调控。

本发明的技术方案是利用改进ABC算法确定热学参数的大坝温度预测方法，所述改进ABC算法为基于交叉运算的全局人工蜂群算法。

所述大坝温度预测方法包括以下步骤：

步骤1：采集大坝不同部位的混凝土温度监测数据；

步骤2：建立各个混凝土浇筑仓的有限元模型，并设定初始条件和边界条件；

步骤3：确定用于大坝温度预测的热学参数及其取值范围，并建立热学参数的目标函数；

步骤4：在热学参数的取值范围内随机产生一组初始热学参数值并代入有限元模型进行温度仿真模拟计算，得到混凝土浇筑仓内节点处的温度值，与节点处混凝土的实际温度一并代入目标函数中，计算得到目标函数值；

步骤5：利用基于交叉运算的全局人工蜂群算法对热学参数进行迭代反演，并利用目标函数判断得到的热学参数值的优劣；

步骤6：重复步骤5，得到满足目标函数值要求的热学参数的取值，并将计算得到的温度值与混凝土实际温度值相比较，达到预期目标后停止迭代反演；

步骤7：根据得到的热学参数值，预测大坝混凝土的温度。

步骤1中，在混凝土浇筑仓内埋设光纤，通过光纤测温获取混凝土内部不同深度处真实温度数据。

步骤2包括以下子步骤：

步骤2-1：在CAD软件中建立起浇筑仓的实体模型，并从CAD软件中提取浇筑仓的坐标值，导入ANSYS软件中，建立浇筑仓的有限元模型，并对该有限元模型进行网格划分；同时对混凝土材料、单元属性进行分配赋值；在ANSYS软件中，选择可支持生死单元的solid70作为混凝土单元；

步骤2-2：有限元模型的初始条件包含混凝土的初始温度，混凝土的初始温度根据现场的实测数据获得；

将浇筑仓底面和横缝侧面作为绝热边界，浇筑仓顶面和上下游表面为第三类温度边界条件。

混凝土与空气接触时，假定经过混凝土表面的热流量q与混凝土表面温度T和气温T_a之差成正比，第三类传热边界条件的表达式如下：

式中q为热流量；λ为导热系数；T_a*为综合等效气温；β为混凝土表面放热系数；T_s为混凝土表面日平均温度；n为法向方向单位矢量。

步骤2-3：边界条件主要包括外界环境温度和混凝土的通水冷却；

一期通水冷却阶段通常会根据实际情况采用多档通水，水温和通水流量会随之发生改变。

考虑冷却通水时混凝土温度的计算式如下：

T(t)＝T_ωi+(T_i-T_ωi)φ_i(t)+θ₀ψ_i(t) (2)

式中T(t)表示混凝土平均温度；e为数学中的自然常数；T_ωi为第i档冷却通水水温；T_i为i-1档通水结束且第i档通水开始时的混凝土温度；φ_i(t)为第i档通水时的水冷函数；θ₀为最终绝热温升；ψ_i(t)为第i档通水时的水冷温升函数；p_i为第i档通水时的水冷参数；t_i为流量或水温改变时刻，当流量或水温改变时，t必须从0开始；t为冷却时间；s、m₁、m₂均为待定的系数。

水冷参数的计算式如下：

k＝2.09-1.35ξ+0.320ξ² (7)

ξ＝λL/c_wρ_wq_w (8)

非金属冷却水管的等效导温系数计算式为：

式中D为等效冷却柱体直径；b为等效冷却柱体半径；S₁、S₂分别为水管布置的水平、铅直间距；ρ_w为水的密度，q_w为通水流量，L为水管长度，c_w为水的比热，λ为混凝土导热系数；a为等效导温系数，c为水管外半径，r₀为水管内半径，λ₁为水管导热系数；ξ表示与水比热、水的密度、通水流量、水管长度有关的中间变量。

步骤3包括以下子步骤：

步骤3-1：确定反演的热学参数及取值范围，反演的热学参数包括混凝土表面放热系数β，混凝土温升规律参数n、混凝土最终绝热温升θ₀；

步骤3-2：建立热学参数反演的目标函数，将混凝土热学参数反演问题转变为最优化问题；

目标函数f(x)的表达式如下：

min f(x)x＝(x₁，x₂，x₃，…，x_D)∈[L_d，U_d] (14)

式中x表示该优化问题的可行解，x^*表示该优化问题的最优解，x为D维向量，D表示优化问题的参数数量；T′_ij、T_ij分别表示第j个位置i时刻的混凝土温度计算值、实测值；p表示测点位置总数量，q为总时长；L_d、U_d分别表示解空间的下界、上界。

步骤4包括以下子步骤：

步骤4-1：混凝土温度求解过程是对混凝土的浇筑次序、大气温度、浇筑温度、水化热的仿真模拟；

ANSYS软件不能直接施加混凝土的水化热，为此将水化热转换为生热率，生热率即单位体积混凝土在单位时间内的生热量，可通过对水泥水化热公式求导得到；计算得到生热率后，再对单元施加生热率并计算温度场。

生热率的计算式如下：

式中HGEN为混凝土水化生热率；Q为混凝土中产生的热量；p为混凝土密度；c为混凝土比热；θ(τ)表示混凝土绝热温升；θ′(τ)表示混凝土绝热温升计算式对龄期τ求导；τ表示混凝土龄期。

根据式(15)，可通过给定的绝热温升表达式来确定混凝土的生热率。

步骤4-2：利用ANSYS软件的“单元生成”和“单元杀死”功能实现分仓浇筑混凝土的模拟仿真；分仓浇筑的模拟仿真前，将坝体单元全部“杀死”；仿真温度计算时，根据浇筑进度依次激活坝体单元；

步骤4-3：通过ANSYS软件中参数化程序设计语言(APDL)将不同时间段的混凝土生热率

施加到浇筑仓的温度场上，再对浇筑仓混凝土进行温度场计算，得到混凝土的计算温度T′_ij。

步骤5中，基于交叉运算的全局人工蜂群算法针对人工蜂群算法收敛速度慢的问题，对解空间的搜索方程式进行如下改进，

式中v_id表示第i个粒子的速度；x_id分别表示第i个粒子在搜索空间中的位置；

表示全局最优解向量中的第d个分量；α为[-1，1]中的随机数；β为[0，C]中的随机数，其中C为非负常数；

式(20)通过

项可平衡人工蜂群算法的探索与开发能力，但在一定程度上降低了全局寻优能力。

针对这个问题，引入遗传算法中的交叉运算操作，对人工蜂群算法进一步改进，

式中v_id′表示引入交叉运算后第i个粒子的速度；cr表示交叉算子，一般取值为0.3～0.6；rand为[0，1]的随机数。

当cr取较小值时，改进的人工蜂群算法的开发能力加强，反之，当cr取较大值时，改进的人工蜂群算法算法的探索能力加强；对应不同的优化问题，cr可取不同的值来使算法达到最好的寻优能力，提高了算法对不同优化问题的适应性。

相比现有技术，本发明的有益效果包括：

1)本发明在考虑冷却通水和外部环境气温影响的基础上，根据混凝土仓中已埋设光纤的测温数据，将改进ABC算法应用到混凝土热学参数反演计算，减小了计算的热学参数值与施工期混凝土真实热学参数值的误差。

2)反演结果对于阐明温度变化与热学参数之间的关系具有重要意义，再根据计算得到的热学参数值预测大坝温度，提高了大坝温度预测的精度，有利于对大坝实施精细化的实时温度调控。

3)本发明针对ABC算法开发能力弱的特点，通过在ABC算法解空间的搜寻过程中增加当前全局最优值的引导作用以及在引领蜂领域搜索过程中增加交叉操作，增加了ABC算法的收敛速度，提高了ABC算法的全局寻优能力。

附图说明

下面结合附图和实施例对本发明作进一步说明。

图1为本发明实施例的分布式光纤测温***的示意图。

图2为本发明实施例的热学参数反演框架的示意图。

图3为本发明实施例的混凝土浇筑仓的分布式光纤典型布置图。

图4为本发明实施例的021混凝土浇筑仓的分布式光纤平面布置图。

图5为本发明实施例的混凝土浇筑仓的有限元模型。

图6为本发明实施例的混凝土浇筑仓温度模拟仿真前处理的流程图。

图7为本发明实施例的混凝土浇筑仓通水水温变化曲线图。

图8为本发明实施例的混凝土浇筑仓通水流量变化曲线图。

图9为本发明实施例的混凝土浇筑仓的环境气温变化曲线图。

图10为本发明实施例的ABC算法的目标函数随迭代次数变化曲线图。

图11为本发明实施例的不同蜜源数量的目标函数随迭代次数变化曲线图。

图12为本发明实施例的021混凝土浇筑仓的温度预测值与实测值的对比图。

图13为本发明实施例的022混凝土浇筑仓的温度预测值与实测值的对比图。

图14为本发明实施例的CGABC算法分量交叉运算的示意图。

具体实施方式

实施例以处于四川与云南交界处的白鹤滩混凝土双曲拱坝为例，该大坝是金沙江下游梯级开发中的第二个梯级水电站，共分为31个坝段。该坝地处亚热带季风区，坝址区多年平均气温21.9℃，极端气温温差大、昼夜温差变化明显。其中极端最高气温42.7℃，发生在9月份，极端最低气温2.1℃，发生在12月份。平均水温17.4℃，平均风速1.9m/s。

实施例采用的混凝土标号主要为C30、C35、C40。为分析白鹤滩工程中C40混凝土在受外界浇筑环境影响下的温度变化规律，实施例选取27#坝段020-022浇筑仓，020-022浇筑仓为三个由下到上的连续浇筑仓，利用改进的人工蜂群算法进行混凝土热学参数反演。

如图1所示，利用改进ABC算法确定热学参数的大坝温度预测方法，包括以下步骤：

步骤1：利用大坝各坝段布设的光纤，测取混凝土温度数据；

实施例中，各坝段的光纤选用“L”型布设方案，浇筑仓内光纤布设线路的起点、终点距离横缝面模板3～5m，距离上、下游面模板3m，如图3～4所示。取每一浇筑仓内的顺河向光纤的温度平均值作为中心点实测温度值。同时在典型测点布设温度计，进行相互验证。

步骤2：对混凝土浇筑仓建立有限元模型，并设立初始条件和边界条件；

步骤2-1：在CAD软件中建立起浇筑仓的实体模型，并从该CAD模型中提取坐标值，导入ANSYS软件中，建立起浇筑仓的有限元模型，并对该有限元模型进行网格划分。同时对混凝土材料、单元属性进行分配赋值。ANSYS软件中，选择可支持生死单元的solid70作为混凝土单元。

步骤2-2：初始条件主要指混凝土的初始温度，由现场的实测数据获得。实施例中，将浇筑仓底面和横缝侧面作为绝热边界，浇筑仓顶面和上下游表面为第三类温度边界条件。

混凝土与空气接触时，假定经过混凝土表面的热流量q与混凝土表面温度T和气温T_a之差成正比，为此建立第三类传热边界条件表达式：

式(1)中q为热流量；λ为导热系数；T_a*为综合等效气温；β为混凝土表面放热系数，单位kJ/(m²·h·℃)；T_s为混凝土表面日平均温度；n为法向方向单位矢量。

步骤2-3：外界条件主要包括外界环境温度和混凝土的通水冷却。其中，一期通水冷却阶段通常采用多档通水，水温和通水流量会随之发生改变。

考虑冷却通水时混凝土平均温度的计算公式如下：

T(t)＝T_ωi+(T_i-T_ωi)φ_i(t)+θ₀ψ_i(t) (2)

式中e为数学中的自然常数，T_ωi为第i档通水水温；T_i为i-1档通水结束且第i档通水开始时的混凝土温度；φ_i(t)为第i档通水时的水冷函数；θ₀为最终绝热温升；ψ_i(t)为第i档通水时水冷温升函数；p_i为第i档通水时水冷参数；t为冷却时间；t_i为流量或水温改变时刻，当流量或水温改变时t必须从0开始计时；s、m₁、m₂均为待定的系数。

水冷参数的计算式如下：

k＝2.09-1.35ξ+0.320ξ² (7)

ξ＝λL/c_wρ_wq_w (8)

非金属冷却水管的等效导温系数计算式为：

式(5)～(8)中：k、a′、D为常数；D为等效冷却柱体直径；b为等效冷却柱体半径；S₁、S₂分别为水管布置的水平、铅直间距；ρ_w为水的密度，q_w为通水流量，L为水管长度，c_w为水的比热，λ为混凝土导热系数；a′为等效导温系数，根据已知的等效冷却柱体半径b、水管外半径c、水管内半径r₀、混凝土导热系数λ、水管导热系数λ₁，带入式(9)求得。

实施例中，通过ANSYS有限元软件建立拱坝27#坝段020～022浇筑仓的有限元模型，这三个连续浇筑仓混凝土标号均为C40且级配相同，每个浇筑仓均分为6个0.5m厚的坯层，仓厚为3m。将浇筑仓的有限元模型通过六面体八节点等参单元进行网格划分，总计74976个单元，83655个节点，如图5所示。C40混凝土密度ρ为2601kg/m³，比热c为0.851kJ/(kg℃)，导温系数α为0.0078m²/d。

将浇筑仓横缝侧面作为绝热边界，浇筑仓顶面和上下游表面为第三类温度边界条件，冷却通水的时间、水温、流量采用时间通水方案，外界环境气温采用实测0.5d平均气温，为了消除019浇筑仓对020浇筑仓混凝土温度的影响，020浇筑仓混凝土实测温度不参与反演计算。按照实际施工中的从下往上的浇筑顺序，利用ANSYS软件的热分析过程中的生死单元，即开始计算前“杀死”所有单元，020浇筑仓在浇筑时只激活020浇筑仓单元，021浇筑仓在浇筑时只激活020、021浇筑仓单元，022浇筑仓浇筑时激活3个浇筑仓所有单元。温度荷载加载时间步长取0.5d，总时间为30d。

ANSYS软件的温度仿真前处理流程图如图6所示。水管冷却效果采用等效热传导法模拟计算，实际水管通水的冷却水温、通水流量和外界环境气温信息分别如图7、8、9所示。

步骤3：确定待反演的热学参数及取值范围；

步骤3-1：将混凝土表面放热系数β，混凝土温升规律参数n、混凝土最终绝热温升θ₀作为待反演参数。考虑到水泥的水化热是影响混凝土温度应力的一个重要因素，而实际上温度场计算中用的是混凝土的绝热温升θ，绝热温升的计算需要用到混凝土的热学参数。

其中混凝土绝热温升表达式常用双曲线式，双曲线式的表达式如下：

式中θ(τ)表示混凝土绝热温升函数，τ表示混凝土龄期。

虽然双曲线式与试验资料吻合度高，但不便运算，由此朱伯芳提出了用于计算混凝土绝热温升的组合指数式：

式(11)的参数之间的关系如下：

式中θ₀为混凝土最终绝热温升，τ为混凝土龄期，s、m₁、m₂为3个待定的系数，实施例中s取值为0.60。

实施例中，选取组合指数式作为混凝土绝热温升表达式。进一步，考虑到混凝土主要热学参数有等效表面放热系数β、温升规律参数n、最终绝热温升θ₀、导热系数λ、导温系数α、混凝土的密度ρ和比热c。其中混凝土密度ρ、比热c和导温系数α可通过室内试验获得比较精确的结果；导温系数λ可根据公式α＝λ/cρ来确定。而混凝土浇筑后需表面覆盖保温被和流水养护，且横缝直接暴露在空气中，需要对等效表面放热系数β进行反演。

因此，本发明选择最终绝热温升θ₀、温升规律参数n、等效表面放热系数β这三个热学参数作为待反演参数。

参考该区混凝土设计值以及长江科学院室内实验值，待反演的热学参数的范围：θ₀∈[18，30]、n∈[1.2，6]、β∈[100，1200]。

步骤3-2：建立目标函数；

为了将经人工蜂群算法反演后的热学参数计算出典型节点的温度值与实际得到的该节点温度值的差值达到最小，即当x^＊＝[θ₀ ^*，n^*，β^*]时，f(x^＊)＝minf(x)，建立以下目标函数表达式(13)，从而将混凝土热学参数反演问题转变为求f(x)函数的最优化问题。

min f(x)x＝(x₁，x₂，x₃，…，x_D)∈[L_d，U_d] (14)

式中f(x)为智能优化问题的目标函数；x为该优化问题的可行解，x^＊为该优化问题的最优解，每个解都是一个D维向量；D为智能优化问题的参数个数；T′_ij、T_ij分别代表i时刻第j个位置的混凝土温度计算值、实测值；p为测点位置总数，q为总时长。L_d、U_d分别为解空间的下界、上界。

步骤4：将取值范围内随机选取的一组热学参数值带入ANSYS中进行温度仿真模拟计算，得到混凝土浇筑仓各节点不同时间的温度计算值T′_ij。再将混凝土的实际温度和计算温度带入目标函数(13)进行计算，得到目标函数值。

步骤4-1：混凝土温度求解过程的主要任务是要对混凝土的浇筑次序、大气温度、浇筑温度、水化热等的仿真模拟。但由于ANSYS软件不能直接施加混凝土的水化热，必须将水化热转换为生热率HGEN，即单位体积混凝土在单位时间内的生热量，可通过对水泥水化热公式求导得到。当计算出生热速率后，再对单元施加生热速率并计算温度场。

生热速率的计算公式为：

其中，HGEN为混凝土水化生热率；Q为混凝土中产生的热量；p为混凝土密度，单位Kg/m³；c为混凝土比热，单位kJ/(kg·℃)。

根据式(15)，就可以通过给定的绝热温升表达式来确定混凝土的生热率。

步骤4-2：由于混凝土浇筑仓并非一次性浇筑完成，而是要分仓浇筑。利用ANSYS中的“单元生死”功能来实现这一过程。计算前，将坝体单元全部“杀死”。计算时，再根据浇筑进度依次激活坝体单元即可。

步骤4-3：通过ANSYS软件中参数化程序设计语言(ANSYS Parametric DesignLanguage，APDL)将不同时间段时的混凝土产热速率施加到温度场上，再对浇筑仓混凝土进行温度场计算，得到混凝土的计算温度T′_ij。

步骤5：通过改进的全局人工蜂群算法将需要反演的热学参数进行迭代，重复步骤4的过程，依靠改进的全局人工蜂群算法的寻优能力得到满足目标函数值要求的热学参数的取值。

人工蜂群算法ABC的标准流程如下：

1)初始化阶段；

将待反演参数的取值范围作为解空间，在解空间内随机产生SN个蜜源，每个蜜源x_i(x_i1，x_i2，…，x_iD)^T(i＝1，2，…，SN)的位置即为一个潜在解，每个蜜源同时只安排一只引领蜂采蜜，其中引领蜂数量＝跟随蜂数量＝SN/2。蜜源i＝1，2，…，SN的初始位置按照下式随机产生：

x_id＝L_d+rand(0，1)(U_d-L_d) (16)

式中rand(0，1)表示0～1之间的随机数。

同时初始化最大开采次数，以便于后面跳出局部最优解。

2)引领蜂阶段；

引领蜂根据式(17)在蜜源i的邻域搜索并产生一个新的蜜源位置vi，依据式(18)计算新蜜源的适应度值fitness(v_i)，采用贪婪选择法与fitness(x_i)进行比较，若v_i的函数值优于x_i，则保留新蜜源的位置信息。

式(16)-(18)中，d是[1，D]中的一个随机整数，表示在x_i(x_i1，x_i2，…，x_iD)^T中的随机一维产生扰动；a是加速度系数，通常取值为1；

是[-1，1]中的随机数，代表扰动幅度；j≠i且j∈{1，2，…，SN}，代表在邻域中选出一个不同于蜜源i的蜜源，随着搜索越来越接近最优解，邻域的范围会逐渐减小；x_i代表蜜源i的位置；fitness(x_i)为蜜源i的适应度值；f(x_i)为蜜源i的函数值。

3)跟随蜂阶段

跟随蜂通过式(18)计算跟随概率p_i，并采用轮盘赌方法选择蜜源，即在[0，1]内随机产生一个随机数r，如果r＜p_i，则跟随蜂选择蜜源i，进行与之前引领蜂一样的步骤，即转到式(16)在邻域内产生新蜜源并进行贪婪选择；若r＞p_i，则进行下一步。

4)侦察蜂阶段

判断是否在最大迭代次数内到达目标函数的要求，若达到，则循环结束，输出最优解；否则放弃该蜜源，在该蜜源采蜜的引领蜂角色转变为侦查蜂，重新返回式(16)随机产生一个新蜜源代替该蜜源，循环继续，直到产生最优解。此时的最优解即反演参数的最优值。

实施例的全局人工蜂群算法(Global Artificial Bee Colony algorithm，GABC)如下：

在步骤1)到步骤4)的算法循环过程中，并没有对当前的全局最优值记忆进行引导，只有随机的一个分量x_i，导致ABC算法开发能力不足，存在收敛速度放慢的问题。针对该问题，本发明选择借鉴粒子群算法中的全局最优解(gbest)来改进标准人工蜂群算法中的搜索方程式(17)，

式(20)中

为算法目前全局最优解向量中的第d个分量；α为[-1，1]中的随机数；β为[0，C]中的随机数，其中C是一个非负常数，当C取0时式(20)等于式(17)，对C取不同值时的算法性能分析试验，C取1.5时算法性能较好，算法的全局寻优能力加强。

式(20)相比于式(17)，通过

项可平衡ABC算法的探索与开发能力，但在一定程度上降低了全局寻优能力，针对这个问题，本发明引入遗传算法中的交叉运算操作，对GABC算法进一步改进，提出基于交叉运算的全局人工蜂群算法(Cross-GlobalArtificial Bee Colony algorithm，CGABC)。

交叉运算是将两个父体的部分码值按位互相交换得到新的个体，在算法中首先每个分量都产生一个0～1的均匀随机数rand，如果rand＜设定值cr，则接受目标分量，否则保留当前分量，如图14所示。

本发明将GABC算法与遗传算法中的二项交叉结合，通过让引领蜂邻域搜索后与当前全局最优值进行交叉操作，即在式(20)后面加入式(21)来提高算法全局寻优能力：

式(21)中cr为算法中设定的交叉算子，一般取值为0.3～0.6；rand为[0，1]之间的一个均匀随机数，当cr取较小值时算法的开发能力加强，反之，当cr取较大值时算法的探索能力加强。

对应不同的优化问题cr可以取不同的值来使算法达到最好的寻优能力，提高了算法对不同优化问题的适应性。

图10给出了标准ABC算法与实施例的CGABC算法在C40混凝土热学参数反演计算中相同迭代次数时的目标函数对比值，可看出两种算法随着迭代次数增加目标函数值都逐渐减小，但CGABC算法相比ABC算法收敛速度更快，在相同迭代次数时目标函数值大部分都比ABC算法更小即反演计算结果更优，说明在热学参数反演时CGABC算法性能更好。

图11给出了改进后的CGABC算法在C40混凝土热学参数反演计算中不同群体规模时的混凝土热学参数反演性能，设定最大迭代次数maxCycle＝200，最大开采次数limit＝100，维数D＝3，通过选取不同的人工蜂群蜜源数量SN进行参数反演计算，得到蜜源数量分别为20、40、60时目标函数值随迭代次数的变化。

步骤6：将反演得到的热学参数带回进行反分析；

实施例中，蜜源数量＝40时反演得到的C40混凝土热学参数：θ₀＝24℃，n＝2.86，β＝489.2kJ/(m²·d·℃)，如表1所示，并将其代入27#坝段021、022浇筑仓的温度场进行仿真计算，得到混凝土不同龄期的计算温度值，021浇筑仓、022浇筑仓的计算温度与实测温度变化趋势对比分别如图12、13所示。从图12和图13可看出021、022浇筑仓混凝土温度实测值与计算值吻合较好，变化趋势一致，其中022浇筑仓的个别时间点出现差异是因为实际施工过程中冷却通水流量不平稳。

表1热学参数反演值表

步骤7：根据得到的热学参数值，预测大坝混凝土的温度。

为了验证本发明的CGABC算法相比改进前的ABC和GABC算法具有更强的寻优能力，实施例中，选取Rastrigin、Sphere、Rosenbrock、Ackely、Schaffer五个常用的性能测试函数来对ABC、GABC、CGABC三种算法分别进行寻优性能测试，其中Sphere函数是单模态函数，其余四个测试函数是多模态函数，所有测试函数的理论最优值为0。这五个性能测试函数的表达式和搜索范围如表2所示。

表2性能测试函数表

在性能对比测试中，蜜源数量＝50，引领蜂数量＝跟随蜂数量＝25，最大迭代次数＝5000，最大开采次数＝300，Rosenbrock函数的维数D＝2，Sphere、Rastrigin、Ackely、Schaffer函数的维数D均取值为30，CGABC算法中的交叉算子取0.5时优化性能较好，所有测试函数在上述设定参数下分别用每种算法单独重复运行10次，算法程序均由matlab语言编写，获得每个函数的平均值mean和标准差SD。

标准人工蜂群算法ABC、改进全局人工蜂群算法GABC和基于交叉运算的全局人工蜂群算法CGABC在5个不同的测试函数时的性能对比如表3所示。从表3的测试结果可以看出CGABC算法相比ABC、GABC算法，在设定的相同条件下寻优能力更强，CGABC算法通过在ABC算法的步骤中引入全局最优解gbest，加快了算法的收敛速度，并与交叉操作结合，显著提高了算法的搜索广度与深度。

表3三种人工蜂群算法的对比表

综上，说明人工蜂群算法能够适用于大坝热学参数反演，反演结果具有较高的可靠性，，从而可为后续大坝混凝土温度实时控制奠定基础。

Claims (6)

1.利用改进ABC算法确定热学参数的大坝温度预测方法，其特征在于，包括以下步骤：

步骤1：采集大坝不同部位的混凝土温度监测数据；

步骤2：建立各个混凝土浇筑仓的有限元模型，并设定初始条件和边界条件；

步骤3：确定用于大坝温度预测的热学参数及其取值范围，并建立热学参数的目标函数；

步骤4：在热学参数的取值范围内随机产生一组初始热学参数值并代入有限元模型进行温度仿真模拟计算，得到混凝土浇筑仓内节点处的温度值，与节点处混凝土的实际温度一并代入目标函数中，计算得到目标函数值；

步骤5：利用基于交叉运算的全局人工蜂群算法对热学参数进行迭代反演，并利用目标函数判断得到的热学参数值的优劣；

步骤6：重复步骤5，得到满足目标函数值要求的热学参数的取值，并将计算得到的温度值与混凝土实际温度值相比较，达到预期目标后停止迭代反演；

步骤7：根据得到的热学参数值，预测大坝混凝土的温度。

2.根据权利要求1所述的大坝温度预测方法，其特征在于，步骤1中，在混凝土浇筑仓内埋设光纤，通过光纤测温获取混凝土内部不同深度处真实温度数据。

3.根据权利要求1所述的大坝温度预测方法，其特征在于，步骤2包括以下子步骤：

步骤2-1：在CAD软件中建立起浇筑仓的实体模型，并从CAD软件中提取浇筑仓的坐标值，导入ANSYS软件中，建立浇筑仓的有限元模型，并对该有限元模型进行网格划分；同时对混凝土材料、单元属性进行分配赋值；

步骤2-2：有限元模型的初始条件包含混凝土的初始温度，混凝土的初始温度根据现场的实测数据获得；

将浇筑仓底面和横缝侧面作为绝热边界，浇筑仓顶面和上下游表面为第三类温度边界条件；

混凝土与空气接触时，假定经过混凝土表面的热流量q与混凝土表面温度T和气温T_a之差成正比，第三类传热边界条件的表达式如下：

式中q为热流量；λ为导热系数；

为综合等效气温；β为混凝土表面放热系数；T_s为混凝土表面日平均温度；n为法向方向单位矢量；

步骤2-3：边界条件主要包括外界环境温度和混凝土的通水冷却；

考虑冷却通水时混凝土温度的计算式如下：

T(t)＝T_ωi+(T_i-T_ωi)φ_i(t)+θ₀ψ_i(t) (2)

式中T(t)表示混凝土平均温度；e为自然常数；T_ωi为第i档冷却通水水温；T_i为i-1档通水结束且第i档通水开始时的混凝土温度；φ_i(t)为第i档通水时的水冷函数；θ₀为最终绝热温升；ψ_i(t)为第i档通水时的水冷温升函数；p_i为第i档通水时的水冷参数；t_i为流量或水温改变时刻；t为冷却时间；s、m₁、m₂均为待定的系数；

水冷参数的计算式如下：

k＝2.09-1.35ξ+0.320ξ² (7)

ξ＝λL/c_wρ_wq_w (8)

非金属冷却水管的等效导温系数计算式为：

式中D为等效冷却柱体直径；b为等效冷却柱体半径；S₁、S₂分别为水管布置的水平、铅直间距；ρ_w为水的密度，q_w为通水流量，L为水管长度，c_w为水的比热，λ为混凝土导热系数；a'为等效导温系数，c为水管外半径，r₀为水管内半径，λ₁为水管导热系数；ξ表示与水比热、水的密度、通水流量、水管长度有关的中间变量。

4.根据权利要求3所述的大坝温度预测方法，其特征在于，步骤3包括以下子步骤：

步骤3-1：确定反演的热学参数及取值范围，反演的热学参数包括混凝土表面放热系数β，混凝土温升规律参数n、混凝土最终绝热温升θ₀；

步骤3-2：建立热学参数反演的目标函数，将混凝土热学参数反演问题转变为最优化问题；

目标函数f(x)的表达式如下：

minf(x)x＝(x₁,x₂,x₃,…,x_D)∈[L_d,U_d] (14)

式中x表示该优化问题的可行解，x^*表示该优化问题的最优解，x为D维向量，D表示优化问题的参数数量；T′_ij、T_ij分别表示第j个位置i时刻的混凝土温度计算值、实测值；p表示测点位置总数量，q为总时长；L_d、U_d分别表示解空间的下界、上界。

5.根据权利要求4所述的大坝温度预测方法，其特征在于，步骤4包括以下子步骤：

步骤4-1：混凝土温度求解过程是对混凝土的浇筑次序、大气温度、浇筑温度、水化热的仿真模拟；

ANSYS软件不能直接施加混凝土的水化热，为此将水化热转换为生热率，生热率即单位体积混凝土在单位时间内的生热量，可通过对水泥水化热公式求导得到；计算得到生热率后，再对单元施加生热率并计算温度场；

生热率的计算式如下：

式中HGEN为混凝土水化生热率；Q为混凝土中产生的热量；p为混凝土密度；c为混凝土比热；θ(τ)表示混凝土绝热温升；θ′(τ)表示混凝土绝热温升函数的导数；τ表示混凝土龄期；

根据式(15)，可通过绝热温升表达式来确定混凝土的生热率；

步骤4-2：利用ANSYS软件的“单元生成”和“单元杀死”功能实现分仓浇筑混凝土的模拟仿真；分仓浇筑的模拟仿真前，将坝体单元全部“杀死”；仿真时，根据浇筑进度依次激活坝体单元；

步骤4-3：通过ANSYS软件中参数化程序设计语言将不同时间段的混凝土生热率施加到浇筑仓的温度场上，再对浇筑仓混凝土进行温度场计算，得到混凝土的计算温度T′_ij。

6.根据权利要求5所述的大坝温度预测方法，其特征在于，步骤5中，基于交叉运算的全局人工蜂群算法针对人工蜂群算法收敛速度慢的问题，对解空间的搜索方程式进行如下改进，

式中v_id表示第i个粒子的速度；x_id表示第i个粒子在搜索空间中的位置；

表示全局最优解向量中的第d个分量；α为[-1,1]中的随机数；β为[0,C]中的随机数，其中C为非负常数；

式(20)通过

项可平衡人工蜂群算法的探索与开发能力，但降低了全局寻优能力，针对这个问题，引入遗传算法中的交叉运算操作，对人工蜂群算法进一步改进，

式中v_id'表示引入交叉运算后第i个粒子的速度；cr表示交叉算子；rand为[0,1]的随机数；

当cr取较小值时，改进的人工蜂群算法的开发能力加强，反之，当cr取较大值时，改进的人工蜂群算法算法的探索能力加强；对应不同的优化问题，cr可取不同的值来使算法达到最好的寻优能力，提高了算法对不同优化问题的适应性。

Images