CN110222389B - 基于虚拟厚度的反击式转轮全局特性数值计算方法和系统 - Google Patents

Abstract

本发明公开了一种基于虚拟厚度的反击式转轮全局特性数值计算方法和系统。步骤包括：数据输入；坐标变换；过水断面形成线构建；轴面速度计算；轴面流线计算；S1、S2流面构建；空间相对流线构建及分割；空间相对水流角计算；空间相对流线上各点轴面速度计算；增厚系数试算；虚拟厚度叶片排挤系数计算；欧拉能量数值求解；效率数值求解；CFD数值模拟结果验证；判断精度；增厚系数插值；将计算结果输入到Excel表格进行显示和存储；计算结束。所述系统用于实现对应功能。该数值计算方法和对应的系统提高了计算过程精简性和高效性；采用模块化设计思路，各子模块既可独立运行又相互之间紧密联系，计算流程清晰；适用于工作在非压缩流体中的反击式叶轮旋转机械，具有广泛应用对象，涉及众多行业领域。

Description

基于虚拟厚度的反击式转轮全局特性数值计算方法和系统

技术领域

本发明涉及反击式旋转机械工作特性数值计算领域，尤其涉及一种反击式转轮全局工况下效率和水力特性的快速预测方法。

背景技术

水力发电因其输出电能稳定、调节反应迅速已成为电网保证高质量电能不可或缺的能源结构，因此水电为了高效发挥其调节性能，要求其能量转换的核心元件——转轮工作在偏离最优工况。我国，始于上世纪70年代对转轮性能的试验研究和始于90年代的基于有限元软件的数值模拟研究经过多年发展日趋成熟，但由于与这两种研究方法相关的辅助设备未能得到充分发展，导致的高成本投入、长时间耗资，比如模型试验的长周期、数值模拟的高计算性能要求等，已成为新建水力发电项目初期阶段不可忽视的成本投入。

发明内容

本发明所要解决的技术问题是：提出基于虚拟厚度的反击式转轮全局特性数值计算方法和系统，缩短研发周期，降低成本。

本发明为解决上述技术问题采用以下技术方案：

本发明提出一种基于虚拟厚度的反击式转轮全局特性数值计算方法，进一步地：步骤包括：

S001、数据输入：读取反击式转轮的全局参数和木模图参数；

S002、坐标变换：将笛卡尔坐标系内的空间叶片投影到Z-r平面上，且遵循以下旋转投影公式：

r_i ²＝x_i ²+y_i ²，z_i＝z_i

其中，(x_i，y_i，z_i)为点在笛卡尔坐标系内的坐标，(z_i，r_i)为对应的点在Z-r平面内的坐标；

将笛卡尔坐标系内坐标(x,y,z)变换到S1流面坐标系内坐标(m,θ)，其中m为S1流面的母线，θ为S1流面圆周方位角；流线的坐标变换并无量纲化处理遵循以下公式：

式中：x₀，y₀为流线起始端点坐标分量，(m_i，θ_i)为笛卡尔坐标系内点(x_i，y_i，z_i)在S1流面坐标系内的坐标；

S003、过水断面形成线构建：在反击式转轮轴面流道内作内切于上冠和下环流线的公切圆，圆心为O，切点为A和B，连接OA、OB和AB，作通过A、B并分别与OA、OB相切的圆弧即为过水断面形成线；

S004、轴面速度计算：求得过水断面形成线形成过水断面时的旋转半径，进一步求得对应的过水断面面积和轴面速度；针对多条过水断面形成线，求得轴面速度在轴面流道中的离散分布；所述过水断面是由过水断面形成线绕转轮转轴旋转一周形成的回转面；

S005、轴面流线计算：读取一条流线对应圆弧的圆心坐标、半径、所述圆弧两端点坐标；列出弧上任意一点的非线性方程组；求得弧上该点坐标；针对多条过水断面形成线，求得一条流线上的多个离散点，顺水流方向用样条曲线依次连接所述离散点，得到所述轴面流线；

S006、S1流面构建：将上一步求得的流线绕转轮旋转轴旋转一周，得到对应的S1流面；

S007、S2流面构建：将木模图参数导入三维建模软件Ansys-Gambit，建立由叶片各截面骨线形成的叶片空间骨面，所述骨面是两相邻叶片之间的一个平均流面，即S2流面；

S008、空间相对流线构建及分割：所述空间相对流线为S1流面与S2流面的交线；将该空间流线等分成n份，输出等分节点及端点直角坐标值，共n+1个点，并利用步骤S002转换到S1流面坐标系中；

S009、空间相对水流角计算：即对步骤S008中每一条空间相对流线上n+1个点处的相对水流角进行计算，对各点处流线切线与水平线的夹角进行组合插值，得到空间相对流线上各点处的空间相对水流角β_i；

S010、空间相对流线上各点轴面速度计算：线性插值轴面流线上离散分布的轴面速度，得到空间相对流线上各等分节点的轴面速度；

S011、叶片厚度计算：根据步骤S001读取的木模图参数，以其中的内切圆圆心坐标和直径数据作为插值参考点，对步骤S008中每一条空间相对流线上n+1个节点处叶片在圆周方向的厚度t_i进行高维线性插值计算；

S012、叶片虚拟厚度计算：根据步骤S011，将空间相对流线每个节点处的厚度乘以增厚系数ψ，得节点处的虚拟厚度，如下式所示：

T_i＝ψ_i·t_i

式中：ψ_i表示增厚系数；T_i表示叶片增厚后的虚拟厚度，单位是m；

S013、叶片排挤系数计算：叶片排挤系数的计算通过下式求得：

式中，

表示节点i处叶片排挤系数；N表示转轮叶片个数；

S014、欧拉能量数值求解：

E_Eui＝U_i·V_Ui＝r_iω·(r_iω-W_Ui)＝r_iω·(r_iω-k_iV_micotβ_i)

式中，E_Eui表示节点i处欧拉能量，单位是m²/s²；U_i表示节点i处圆周速度，单位是m/s；V_Ui表示节点i处绝对速度圆周分量，单位是m/s；ω表示转轮旋转角速度，单位是rad/s；W_Ui表示节点i处相对速度圆周分量，单位是m/s；V_mi表示节点i处轴面速度，单位是m/s；

S015、转轮效率解析求解：

效率的求解基于水轮机基本能量方程：

式中，γ表示流体比重，单位是N/m³；H表示转轮进出口净水头，单位是m；η_r1表示转轮效率解析解；Q表示活动导叶开度α下的流量，单位是m³/s，g表示重力加速度，单位是m/s²；U₁表示节点i＝1处圆周速度，单位是m/s；V_U1表示i＝1处绝对速度圆周分量，单位是m/s；U_n+1表示节点i＝n+1处圆周速度，单位是m/s；V_Un+1表示i＝n+1处绝对速度圆周分量，单位是m/s；

S016、CFD数值模拟结果验证：利用软件CFX对工作在与步骤S015相同工况下的转轮进行数值模拟；通过数值模拟得到转轮中流体流动特性以及转轮转轴上的扭矩，根据以下方程求得转轮相应效率：

Mω＝γQHη_r2

式中：M表示转轮转轴扭矩，单位是N·m；

S017、判断精度：对步骤S015和S016所得的转轮效率求解相对误差，若相对误差超出判定标准值，则返回步骤S012，改变增厚系数取值，再依次执行步骤S012～S017；若相对误差小于判定标准值，则执行步骤S018；所述判定标准值预先给定；

S018、令转轮的活动导叶在活动范围内等间隔地取导叶开度的值，每一个导叶开度对应一个流量工况，对所述流量工况按照对应导叶开度由小到大依次编号；间隔地取所述流量工况，依次执行步骤S004～S017，得到所述流量工况下转轮效率的解析解与对应的增厚系数；

S019、增厚系数插值：利用现有的导叶开度与增厚系数数据，对增厚系数ψ进行插值，求得其余导叶开度对应的流量工况下的增厚系数；

S020、根据步骤S019计算所得的增厚系数，依次执行步骤S004～S016，然后计算步骤S015所得转轮效率解析解和步骤S016所得转轮效率模拟解的相对误差；

S021、将计算结果输入到Excel表格进行显示和存储，计算结束。

如前所述的一种基于虚拟厚度的反击式转轮全局特性数值计算方法，进一步地：所述反击式转轮全局参数包括：所述反击式转轮叶片个数、设计水头、额定转速、不同导叶开度及与之对应的流量、流体密度和环境温度；所述反击式转轮的木模图参数包括：所述反击式转轮木模图上各叶片水平截面内切圆圆心坐标、圆与圆切点坐标和内切圆直径。

如前所述的一种基于虚拟厚度的反击式转轮全局特性数值计算方法，进一步地：步骤S005所述的非线性方程组为：

(z_F-z_O′)²+(y_F-y_O′)²＝R²

(z_F-z_A)²+(y_F-y_A)²＝R²+R²-2R²cosα₁

(z_F-z_B)²+(y_F-y_B)²＝R²+R²-2R²cosα₂

其中，所述流线对应圆弧为

其圆心为O′(z_O′，y_O′)，半径为R，

端点A和B的坐标分别为(z_A，y_A)、(z_B，y_B)；F(z_F，y_F)为弧

上任一点；α₁为∠AO′F弧度值，单位是rad；α₂为∠BO′F弧度值，单位是rad。

如前所述的一种基于虚拟厚度的反击式转轮全局特性数值计算方法，进一步地：步骤S009所述求得空间相对水流角的值的过程包括：

根据最优工况下流体相对流动方向与空间相对流线切线方向一致原则，对各点处流线切线与水平线的夹角进行组合插值，求得空间相对流线上各点处的空间相对水流角的值，根据点所在位置不同，分以下三种情况：

i＝1

2≤i≤n

i＝n+1

式中：θ为空间相对流线上各点在S1流面坐标系内的θ轴分量值；β_i为节点i处相对水流角，单位是rad。

如前所述的一种基于虚拟厚度的反击式转轮全局特性数值计算方法，进一步地：步骤S012中只计算流线第一个节点和最后一个节点处的虚拟厚度。

本发明还提出一种基于虚拟厚度的反击式转轮全局特性数值计算系统，进一步地：包括数据采集模块、前处理模块、数值计算模块和显示模块，所述数据采集模块的输出传至前处理模块，所述前处理模块的输出分别传至数值计算模块和显示模块；所述数据采集模块用于采集所述反击式转轮的全局参数和木模图参数；所述前处理模块用于变换坐标、计算轴面流线、计算轴面速度、计算S1流面、S2流面、计算空间相对流线、计算空间相对流线上节点处轴面速度、计算空间相对水流角、试算叶片增厚系数、计算虚拟厚度叶片的排挤系数和增厚系数插值；所述数值计算模块用于由准备参数求得目标参数的计算值，由数值模拟求得目标参数的CFD模拟值，以及求得所述目标参数的计算值和CFD模拟值之间的相对误差；所述显示模块用于存储和显示所述目标参数的计算值及对应的增厚系数、CFD模拟值和相对误差；所述准备参数包括：所述空间相对流线上各节点坐标、所述空间相对流线上各节点处轴面速度、所述空间相对流线上各节点处相对水流角和所述空间相对流线上各节点处虚拟厚度叶片的排挤系数；所述目标参数包括反击式转轮进出口的欧拉能量和水力效率。

如前所述的一种基于虚拟厚度的反击式转轮全局特性数值计算系统，进一步地：

所述数据采集模块包括已知数据库、用户自采集子模块和实时数据库，所述已知数据库用于存储所述反击式转轮的全局参数，所述用户自采集子模块用于采集所述反击式转轮的木模图参数，所述实时数据库用于存储所述已知数据库和所述用户自采集子模块中的数据；

所述前处理模块包括坐标变换子模块、轴面流线计算子模块、轴面速度计算子模块、S1流面计算子模块、S2流面计算子模块、空间相对流线计算子模块、空间相对流线轴面速度计算子模块、空间相对水流角计算子模块、增厚系数试算子模块、虚拟厚度叶片排挤系数计算子模块和增厚系数插值子模块；所述前处理模块中的坐标变换子模块从所述数据采集模块获取所述全局参数和木模图参数；所述坐标变换子模块输出所述全局参数、木模图参数的轴面投影形式和木模图参数在S1流面坐标系下的表现形式，分别传输至轴面流线计算子模块、轴面速度计算子模块、S1流面计算子模块和S2流面计算子模块；所述轴面流线计算子模块、S1流面计算子模块和S2流面计算子模块输出各自模块的计算结果，都传至空间相对流线计算子模块；所述空间相对流线计算子模块的第一路输出连同所述轴面速度计算子模块的结果一起传输至空间相对流线轴面速度计算子模块，第二路输出传至空间相对水流角计算子模块，第三路输出连同增厚系数试算子模块的结果一起传输至虚拟厚度叶片排挤系数计算子模块；增厚系数插值子模块从数值计算模块获取初始输入值，对已求得的导叶开度—增厚系数进行插值，计算结果返回到虚拟厚度叶片排挤系数计算子模块；所述前处理模块用于计算得到求解目标参数所需要的准备参数，并将所述准备参数传输至所述显示模块的数据存储子模块；

所述数值计算模块包括数值计算子模块和计算结果验证子模块，其中，所述数值计算子模块的输入是来自空间相对流线计算子模块、空间相对流线轴面速度计算子模块、空间相对水流角计算子模块和虚拟厚度叶片排挤系数计算子模块的准备参数；所述数值计算子模块的输出为目标参数的计算值和增厚系数；所述目标参数的计算值输出至计算结果验证子模块；所述计算结果验证子模块通过CFD数值模拟求得所述目标参数的模拟值，并且求得所述目标参数的计算值和模拟值之间的相对误差；若所述相对误差超出判定标准值，则改变增厚系数试算子模块中增厚系数的初始赋值，返回到虚拟厚度叶片排挤系数计算子模块并再次计算；当相对误差小于判定标准值，即认为所述目标参数的计算值是最终解析解，且所述增厚系数满足误差精度要求；选取一定数量的离散工况，每隔一个工况执行上述过程，将每个工况的增厚系数和对应的导叶开度传递至所述前处理模块的增厚系数插值子模块，插值得到剩余工况对应的增厚系数，再依次执行上述子模块，得剩余工况下所述目标参数计算值和相对误差；将所述计算值、增厚系数、CFD模拟值和相对误差一路传输至所述显示模块的计算结果显示子模块进行显示和输出，另一路传输至所述显示模块的数据存储子模块进行存储；

所述显示模块包括计算结果显示子模块和数据存储子模块，所述计算结果显示子模块用于显示和输出所述目标参数的最终计算值、所述增厚系数、所述CFD模拟值和所述相对误差，所述数据存储子模块用于存储所述准备参数、所述目标参数计算值、所述增厚系数、所述CFD模拟值和所述相对误差。

如前所述的一种基于虚拟厚度的反击式转轮全局特性数值计算系统，进一步地：所述空间相对流线上各节点处相对水流角计算基于流体相对流动方向与所述空间相对流线切向方向一致的假设，计算结果为所述空间相对流线与水平线之间夹角的组合插值。

如前所述的一种基于虚拟厚度的反击式转轮全局特性数值计算系统，进一步地：所述空间相对流线上各点节处叶片厚度根据所述反击式转轮叶片个数、所述内切圆直径通过高维线性插值求得；所述空间相对流线上各节点处叶片厚度排挤系数的计算依据流道输送能力不变原则。

本发明采用以上技术方案与现有技术相比，具有以下技术效果：

该方法将解析法和数值模拟法结合，在部分工况的解析解求解过程中，参数取值试算的步骤通过插值，使参数的初始取值更精确，从而使达到要求的试算次数大大减少。

该数值解析方法一定程度上避免了基于有限元方法的大量三维数值模拟计算，计算高效，时间投入大大缩短；计算方法输入数据简单且容易获得，提高了计算过程精简性和快速性；前处理模块基于Excel、AutoCAD、Matlab和Ansys-Gambit三维建模软件，通过简单指令完成所述全局参数、所述木模图参数和所述准备参数的获取，处理流程易于掌握，针对不同用户适应性强；计算方法采用模块化设计思路，各子模块既可独立运行又相互之间紧密联系，计算流程清晰；计算方法适用于工作在非压缩流体中的反击式叶轮旋转机械，具有广泛应用对象，涉及众多行业领域。

附图说明

图1是本发明一个实施例的数值计算系统结构示意图；

图2是本发明一个实施例的数值计算方法步骤流程图；

图3是本发明一个实施例的转轮叶片空间骨面(S2流面)和轴面投影示意图；

图4是本发明一个实施例的S1流面坐标系示意图；

图5是本发明一个实施例的轴面流道某一过水断面形成线示意图；

图6是本发明一个实施例的轴面流线上某一点坐标计算示意图；

图7是本发明一个实施例的空间相对流线示意图；

图8是本发明一个实施例的空间相对流线相对水流角计算示意图；

图9是本发明一个实施例的叶片有限个数和叶片厚度影响示意图；

图10是本发明一个实施例的叶片虚拟厚度影响示意图。

具体实施方式

下面结合附图对本发明的技术方案做进一步的详细说明：

本技术领域技术人员可以理解的是，除非另外定义，这里使用的所有术语(包括技术术语和科学术语)具有与本发明所属领域中的普通技术人员的一般理解相同的意义。还应该理解的是，诸如通用字典中定义的那些术语应该被理解为具有与现有技术的上下文中的意义一致的意义，并且除非像这里一样定义，不会用理想化或过于正式的含义来解释。

如图1所示，本发明的目的在于提供一种基于虚拟厚度的反击式转轮在其运行的全局范围内进出口欧拉能量和效率的数值计算方法，该方法包括以下模块：

数据采集模块，包括已知数据库、用户自采集子模块以及用以临时存储已知数据库和用户自采集子模块中数据的实时数据库；已知数据库用Excel表格存储反击式转轮全局参数，包括转轮叶片个数、设计水头、额定转速、不同导叶开度及与之对应的流量、流体密度和环境温度；用户自采集子模块用Excel表格存储反击式转轮木模图上相关参数，包括转轮木模图上各叶片水平截面内切圆圆心坐标、圆与圆切点坐标和内切圆直径；这些数据将临时存储到实时数据库的Excel表格中，等待前处理模块的调用。

前处理模块，包括坐标变换子模块、轴面流线计算子模块、轴面速度计算子模块、S1流面计算子模块、S2流面计算子模块、空间相对流线计算子模块、空间相对流线轴面速度计算子模块、空间相对水流角计算子模块、增厚系数试算子模块、虚拟厚度叶片排挤系数计算子模块和增厚系数插值子模块；一方面前处理模块从数据采集模块调用数据，经过一系列前处理子模块得到目标参数求解所需要的准备参数，包括空间相对流线上各节点坐标、各节点处轴面速度、各节点处相对水流角、各节点处虚拟厚度叶片的排挤系数；另一方面，前处理模块担负增厚系数求解和插值的功能；这些数据将存储到数据存储子模块的Excel表格中，等待数值计算模块的调用。

数值计算模块，包括欧拉能量及效率计算子模块和计算结果验证子模块；数值计算模块从数据存储子模块调用准备参数，根据旋转机械流动速度三角形求解空间相对流线上各点欧拉能量，并建立水流势能—机械能转换平衡关系式，求解特定工况下转轮的效率；计算结果验证子模块利用计算流体动力学(CFD)数值模拟结果对实施本发明所得目标参数计算值进行验证，通过调整增厚系数试算子模块中增厚系数值，使得两者之间相对误差小于判定标准值。选取一定数量的离散工况，及相应的导叶开度和流量，每隔一个工况执行上述过程，将每个工况的增厚系数和对应的导叶开度传输到前处理模块的增厚系数插值子模块，插值得到未执行上述过程工况对应的增厚系数值，再依次执行上述步骤。最终，将所有工况目标参数计算值、增厚系数、CFD模拟值和相对误差值输入到显示模块的计算结果显示子模块的Excel表格中，同时存储到数据存储子模块的Excel表格中，目标参数包括转轮进出口欧拉能量和效率。

显示模块，包括计算结果显示子模块和数据存储子模块，计算结果显示子模块用于显示和输出目标参数计算值、增厚系数、CFD模拟值和相对误差值，数据存储子模块用于存储准备参数、目标参数计算值、增厚系数、CFD模拟值和相对误差值。

图2是实施本发明计算方法采用技术方案的具体流程图，包括以下步骤：

S001、数据输入，读取全局参数和相关参数；

S002、坐标变换，轴面投影。如图3所示，将空间叶片投影到Z-r/Y平面上，遵循以下旋转投影公式：

r_i ²＝x_i ²+y_i ²，z_i＝z_i

S003、坐标变换，笛卡尔坐标系内坐标(x,y,z)变换到S1流面坐标系内坐标(m,θ)。S1流面坐标系如图4所示，m为S1流面的母线，θ为S1流面圆周方位角，坐标变换并无量纲化处理遵循以下公式：

式中：x₀，y₀——流线起始端点坐标分量；

S004、过水断面形成线构建，如图5所示，在反击式转轮轴面流道内作内切于上冠和下环流线的公切圆，圆心为O，切点为A和B，链接OA、OB和AB，作通过A、B并分别与OA、OB相切的圆弧即为过水断面形成线，其长度可由经验公式计算得：

式中：

——过水断面形成线弧AB长，m；

——弦AB长，m；

——内切圆半径m；

S005、轴面速度计算，基于转轮一元理论设计方法，即轴面速度在同一过水断面上均匀分布。过水断面是由过水断面形成线绕转轮转轴旋转一周形成的回转面。如图5所示，过圆心O作弦

的垂线垂足为点C，将OC三等分，取一点D，使得：

点D到旋转轴距离

即为过水断面形成线

形成过水断面时的旋转半径，则过水断面面积和某一工况下轴面速度大小为：

V_m＝Q/S

式中：

——线段DE长，可在AutoCAD中测量，m；V_′——轴面速度，m/s；Q——特定工况下的流量，m³/s。

针对多条过水断面形成线，运用上述计算流程，可求得轴面速度在轴面流道中的离散分布。选取不同的流量即可求得不同工况下轴面速度在轴面流道中的分布。由于一元理论的设计方法假设转轮叶片无限多、无限薄，此处V_m是在该假设下计算所得；同时，当转轮偏离最优工况时，V_m也会产生相应变化，叶片有限个数、厚度和工况偏离的影响将在步骤S014-S018中叙述；

S006、轴面流线计算，轴面流线计算依据两两流线之间单位时间输送相等流量的原则。如图6所示的轴面流道，以中间流线SL-0.5(上冠流线：SL-0，下环流线：SL-1)的求解为例，圆弧

圆心为O′(z_O′，y_O′)，半径为R，

端点A和B的坐标分别为(z_A，y_A)、(z_B，y_B)，上述值均可在AutoCAD中测量获取。在弧

上任取一点F(z_F，y_F)，设

段和

段输送能力相同，又由于V_m在过水断面上均匀分布，则可列出以下非线性方程组：

(z_F-z_O′)²+(y_F-y_O′)²＝R²

(z_F-z_A)²+(y_F-y_A)²＝R²+R²-2R²cosα₁

(z_F-z_B)²+(y_F-y_B)²＝R²+R²-2R²cosα₂

式中：α₁——∠AO′F弧度值，rad；α₂——∠BO′F弧度值，rad；

S007、步骤S006中非线性方程组可在Matlab中采用fsolve函数求解，最终求得F点坐标，该点既在过水断面形成线

上也在中间流线SL-0.5上。在其他过水断面形成线上运用该方法，可求得流线SL-0.5上的一系列离散点，顺着水流方向用样条曲线依次连接这些点即为SL-0.5流线；类似的，在SL-0、SL-0.5之间可求出流线SL-0.25，在SL-0.5、SL-1之间可求出流线SL-0.75，以此类推可求得流线SL-0.0625，SL-0.125，SL-0.375，SL-0.625和SL-0.875；

S008、S1流面构建，将步骤S007中求得的流线绕转轮旋转轴旋转一周，得到对应的S1流面；

S009、S2流面构建，将木模图参数导入三维建模软件Ansys-Gambit，建立由叶片各截面骨线形成的叶片空间骨面，此骨面可认为是两相邻叶片之间的一个平均流面，即S2流面；

S010、空间相对流线构建及分割，如图7所示，S1流面与S2流面的交线即为空间相对流线，记为SLS-No.。在三维AutoCAD中将该空间流线等分成n份，输出等分节点及端点的笛卡尔坐标值，共n+1个点，并利用步骤S003转换到m-θ坐标系中；

S011、空间相对水流角计算，即对步骤S010中n+1个节点处的相对水流角进行计算。如图8所示，根据叶片无限多、无限薄假设下流体相对流动方向与空间相对流线切线方向一致原则，对各点处流线切线与水平线的夹角进行一定的组合插值，求得空间相对流线上各点处的空间相对水流角的值，根据点所在位置不同，有以下三种情况：

i＝1

2≤i≤n

i＝n+1

式中：β_i——节点i处相对水流角，rad。

S012、空间相对流线上各节点轴面速度计算，在步骤S005中已计算出轴面速度在各个过水断面形成线上的值，同时步骤S006计算出轴面流线与过水断面形成线的交点值，因此得到轴面速度在轴面流线上的离散分布。步骤S010得到空间相对流线上n+1个等分节点，这些节点的m轴坐标分量与空间相对流线对应的轴面流线一致，解释为，若将空间相对流线旋转投影到轴面上，空间相对流线会与其对应的轴面流线重合。在上述“一致性”的原则下，线性插值轴面流线上离散分布的轴面速度，可得到空间相对流线上各等分节点的轴面速度；

S013、步骤S012中的线性插值可在Matlab中采用interp1函数实现；

S014、叶片排挤系数1/k计算，如图9所示，叶片有限个数和厚度对转轮中的流体存在一定的排挤作用，影响速度分量的大小和方向，此处仅考虑对速度大小的影响。根据步骤S001读取的相关参数，以内切圆圆心坐标和直径数据作为插值参考点，对步骤S010中每一条空间相对流线上n+1个节点处叶片在圆周方向的厚度t_i进行高维线性插值计算；

S015、步骤S014中的高维线性插值可在Matlab中采用griddatan函数实现；

S016、叶片虚拟厚度计算，当转轮工作在偏离最优工况时，转轮流道内的速度有相应程度的改变，如图10中虚线所示，此处通过改变叶片厚度来反映此影响，定量地提出叶片虚拟厚度概念，即为叶片厚度乘以一个增厚系数ψ，如下式所示：

T_i＝ψ_i·t_i

式中：ψ_i表示节点i处增厚系数；T_i表示节点i处叶片增厚后的虚拟厚度，单位是m；

S017、步骤S016中对空间流线进出口节点处的厚度进行增厚处理即可；

S018、叶片排挤系数的计算遵循考虑叶片排挤影响前后，转轮流道输送能力不变的原则，可通过下式求得：

式中：

——节点i处叶片排挤系数；N——转轮叶片个数；

S019：数值求解，欧拉能量求解：

E_Eui＝U_i·V_Ui＝r_iω·(r_iω-W_Ui)＝r_iω·(r_iω-k_iV_micotβ_i)

式中：E_Eui——节点i处欧拉能量，m²/s²；U_i——节点i处圆周速度，m/s；V_Ui——节点i处绝对速度圆周分量，m/s；ω——转轮旋转角速度，rad/s；W_Ui——节点i处相对速度圆周分量，m/s；V_mi——节点i处轴面速度，m/s；

S020、数值求解，最高效率求解基于水轮机基本能量方程：

式中：γ——流体比重，N/m³；H——转轮进出口净水头，m；η_r——转轮效率；

S021、CFD数值模拟结果验证，利用软件CFX在相同工况下的数值模拟结果对实施本发明计算方法所得结果进行验证。通过数值模拟得到转轮中流体流动特性以及转轮转轴上的扭矩，根据以下方程求得转轮相应效率：

Mω＝γQHη_r2

式中：M——转轮转轴扭矩，N·m。

S022、判断本发明数值计算结果是否满足精度要求，对步骤S020和S021所得的转轮效率求解相对误差，若相对误差超出判定标准值，则返回步骤S016，改变增厚系数，依次执行步骤S016～S022；若相对误差小于判定标准值，则执行步骤S023；

S023、选取一定数量的转轮运行工况点，每隔一个工况执行上述步骤，得到增厚系数和对应的导叶开度；

S024、增厚系数插值，对步骤S023所得的导叶开度与增厚系数(α～ψ)数据进行插值，求得其他导叶开度对应的增厚系数值；

S025、步骤S024中插值根据已知参考点个数可选择不同阶数的多项式插值；

S026、对步骤S024所得新的增厚系数和对应的导叶开度依次执行步骤S005～S021，再计算步骤S020和S021所得的转轮效率相对误差；

S027、将计算结果输入到Excel表格进行显示和存储，计算结束。

实施例一

本实施例中，通过上述步骤，计算结果部分数据报告如表1所示，主要包括部分全局参数值、计算效率值、增厚系数初始值和插值、CFD数值模拟效率值和相对误差。

表1

上述技术方案实施步骤，可借助Excel、AutoCAD、Matlab和Ansys三维建模软件Gambit实现所有操作，操作指令简单，对计算硬件要求低且用户适应性强；计算流程经过初次计算编辑后，可形成固定计算模板，便于以后集成工作的开展。

以上所述仅是本发明的部分实施方式，应当指出，对于本技术领域的普通技术人员来说，在不脱离本发明原理的前提下，还可以做出若干改进和润饰，这些改进和润饰也应视为本发明的保护范围。

Claims (7)

1.一种基于虚拟厚度的反击式转轮全局特性数值计算方法，其特征在于：

S001、数据输入：读取反击式转轮的全局参数和木模图参数；

S002、坐标变换：将笛卡尔坐标系内的空间叶片投影到Z-r平面上，且遵循以下旋转投影公式：

r_i ²＝x_i ²+y_i ²，z_i＝z_i

其中，(x_i，y_i，z_i)为点在笛卡尔坐标系内的坐标，(z_i，r_i)为对应的点在Z-r平面内的坐标；

将笛卡尔坐标系内坐标变换到S1流面坐标系内，其中m为S1流面的母线，θ为S1流面圆周方位角；流线的坐标变换并无量纲化处理遵循以下公式：

式中：x₀，y₀为流线起始端点坐标分量，(m_i，θ_i)为笛卡尔坐标系内点(x_i，y_i，z_i)在S1流面坐标系内的坐标；

S003、过水断面形成线构建：在反击式转轮轴面流道内作内切于上冠和下环流线的公切圆，圆心为O，切点为A和B，连接OA、OB和AB，作通过A、B并分别与OA、OB相切的圆弧即为过水断面形成线；

S004、轴面速度计算：求得过水断面形成线形成过水断面时的旋转半径，然后求得对应的过水断面面积和轴面速度；针对多条过水断面形成线，求得轴面速度在轴面流道中的离散分布；所述过水断面是由过水断面形成线绕转轮转轴旋转一周形成的回转面；

S005、轴面流线计算：读取一条流线对应圆弧的圆心坐标、半径、所述圆弧两端点坐标；列出弧上任意一点的非线性方程组；求得弧上该点坐标；针对多条过水断面形成线，求得一条流线上的多个离散点，顺水流方向用样条曲线依次连接所述离散点，得到所述轴面流线；

S006、S1流面构建：将上一步求得的流线绕转轮旋转轴旋转一周，得到对应的S1流面；

S007、S2流面构建：将木模图参数导入三维建模软件Ansys-Gambit，建立由叶片各截面骨线形成的叶片空间骨面，所述骨面是两相邻叶片之间的一个平均流面，即S2流面；

S008、空间相对流线构建及分割：所述空间相对流线为S1流面与S2流面的交线；将该空间流线等分成n份，输出等分节点及端点的直角坐标系的坐标值，共n+1个点，并将所述点的坐标转换到S1流面坐标系中；

S009、空间相对水流角计算：对步骤S008中每一条空间相对流线上所述n+1个点处的相对水流角进行计算，对各点处流线切线与水平线的夹角进行组合插值，得到空间相对流线上各点处的空间相对水流角βi；

S010、空间相对流线上各点轴面速度计算：线性插值轴面流线上离散分布的轴面速度，得到空间相对流线上各等分节点的轴面速度；

S011、叶片厚度计算：根据步骤S001读取的木模图参数，以其中的内切圆圆心坐标和直径数据作为插值参考点，对步骤S008中每一条空间相对流线上n+1个节点处叶片在圆周方向的厚度ti进行高维线性插值计算；

S012、叶片虚拟厚度计算：根据步骤S011，将空间相对流线每个节点处的厚度乘以增厚系数ψ，得所述节点处的虚拟厚度，如下式所示：

T_i＝ψ·t_i

式中：ψ表示增厚系数；T_i表示叶片增厚后的虚拟厚度，单位是m；i为空间相对流线上等距离取得的节点编号，所述节点编号从靠近所述反击式转轮输入口到靠近所述反击式转轮输出口依次增大；

S013、叶片排挤系数计算：叶片排挤系数的计算通过下式求得：

式中，

表示节点i处叶片排挤系数；N表示转轮叶片个数；

S014、欧拉能量数值求解：

E_Eui＝U_i·V_Ui＝r_iω·(r_iω-W_Ui)＝r_iω·(r_iω-k_iV_micotβ_i)

式中，E_Eui表示节点i处欧拉能量，单位是m²/s²；U_i表示节点i处圆周速度，单位是m/s；V_Ui表示节点i处绝对速度圆周分量，单位是m/s；ω表示转轮旋转角速度，单位是rad/s；W_Ui表示节点i处相对速度圆周分量，单位是m/s；V_mi表示节点i处轴面速度，单位是m/s；

S015、转轮效率解析求解：

所述转轮效率的解析求解基于水轮机基本能量方程：

式中，γ表示流体比重，单位是N/m³；H表示转轮进出口净水头，单位是m；η_r1表示转轮效率解析解；Q表示活动导叶开度α下的流量，单位是m³/s，g表示重力加速度，单位是m/s²；U₁表示节点i＝1处圆周速度，单位是m/s；V_U1表示i＝1处绝对速度圆周分量，单位是m/s；U_n+1表示节点i＝n+1处圆周速度，单位是m/s；V_Un+1表示i＝n+1处绝对速度圆周分量，单位是m/s；

S016、CFD数值模拟结果验证：利用软件CFX对工作在与步骤S015相同工况下的转轮进行数值模拟；通过数值模拟得到转轮中流体流动特性以及转轮转轴上的扭矩，根据以下方程求得转轮相应效率：

Mω＝γQHη_r2

式中：M表示转轮转轴扭矩，单位是N·m，η_r2表示转轮效率模拟解；

S017、判断精度：对步骤S015和S016所得的转轮效率解析解和转轮效率模拟解求解相对误差，若相对误差超出判定标准值，则返回步骤S012，改变增厚系数取值，再依次执行步骤S012～S017，最终得到符合要求的转轮效率解析解与对应的增厚系数；若相对误差小于判定标准值，则执行步骤S018；所述判定标准值预先给定；

S018、令转轮的活动导叶在活动范围内等间隔地取导叶开度的值，每一个导叶开度对应一个流量工况，对所述流量工况按照对应导叶开度由小到大依次编号；间隔地取所述流量工况，依次执行步骤S004～S017，得到所述流量工况下转轮效率的解析解与对应的增厚系数；

S019、增厚系数插值：利用现有的导叶开度与增厚系数数据，对增厚系数ψ进行插值，求得其余导叶开度对应的流量工况下的增厚系数；

S020、根据步骤S019计算所得的增厚系数，依次执行步骤S004～S016，然后计算步骤S015所得转轮效率解析解和步骤S016所得转轮效率模拟解的相对误差；

S021、将计算结果输入到Excel表格进行显示和存储，计算结束。

2.如权利要求1所述的一种基于虚拟厚度的反击式转轮全局特性数值计算方法，其特征在于：所述反击式转轮全局参数包括：所述反击式转轮叶片个数、设计水头、额定转速、不同导叶开度及与之对应的流量、流体密度和环境温度；所述反击式转轮的木模图参数包括：所述反击式转轮木模图上各叶片水平截面内切圆圆心坐标、圆与圆切点坐标和内切圆直径。

3.如权利要求1所述的一种基于虚拟厚度的反击式转轮全局特性数值计算方法，其特征在于：步骤S005所述的非线性方程组为：

(z_F-z_o′)²+(y_F-y_o′)²＝R²

(z_F-z_A)²+(y_F-y_A)²＝R²+R²-2R²cosα₁

(z_F-z_B)²+(y_F-y_B)²＝R²+R²-2R²cosα₂

其中，(z_O′，y_O′)为所述流线对应圆弧

的圆心O′的坐标，R为所述流线对应圆弧

的半径，(z_A，y_A)、(z_B，y_B)分别为

端点A和B的坐标；(z_F，y_F)为弧

上任一点F的坐标；α₁为∠AO′F弧度值，单位是rad；α₂为∠BO′F弧度值，单位是rad。

4.如权利要求1所述的一种基于虚拟厚度的反击式转轮全局特性数值计算方法，其特征在于：步骤S009所述求得空间相对水流角的值的过程包括：

根据最优工况下流体相对流动方向与空间相对流线切线方向一致原则，对各点处流线切线与水平线的夹角进行组合插值，求得空间相对流线上各点处的空间相对水流角的值，根据点所在位置不同，分以下三种情况：

i＝1

2≤i≤n

i＝n+1

式中：θ为空间相对流线上各点在S1流面坐标系内的θ轴分量值；β_i为节点i处相对水流角，单位是rad。

5.如权利要求1所述的一种基于虚拟厚度的反击式转轮全局特性数值计算方法，其特征在于：步骤S012中只计算所述空间相对流线第一个节点和最后一个节点处的虚拟厚度。

6.一种基于虚拟厚度的反击式转轮全局特性数值计算系统，其特征在于：包括数据采集模块、前处理模块、数值计算模块和显示模块，所述数据采集模块的输出传至前处理模块，所述前处理模块的输出分别传至数值计算模块和显示模块；所述数据采集模块用于采集所述反击式转轮的全局参数和木模图参数；所述前处理模块用于变换坐标、计算轴面流线、计算轴面速度、计算S1流面、S2流面、计算空间相对流线、计算空间相对流线上节点处轴面速度、计算空间相对水流角、试算叶片增厚系数、计算虚拟厚度叶片的排挤系数和计算增厚系数插值；所述数值计算模块用于由准备参数求得目标参数的计算值，由数值模拟求得目标参数的模拟值，以及求得所述目标参数的计算值和模拟值之间的相对误差；所述显示模块用于存储和显示所述目标参数的计算值及对应的增厚系数、CFD模拟值和相对误差；所述准备参数包括：所述空间相对流线上各节点坐标、所述空间相对流线上各节点处轴面速度、所述空间相对流线上各节点处相对水流角和所述空间相对流线上各节点处虚拟厚度叶片的排挤系数；所述目标参数包括反击式转轮进出口的欧拉能量和水力效率；

所述数据采集模块包括已知数据库、用户自采集子模块和实时数据库，所述已知数据库用于存储所述反击式转轮的全局参数，所述用户自采集子模块用于采集所述反击式转轮的木模图参数，所述实时数据库用于存储所述已知数据库和所述用户自采集子模块中的数据；

所述前处理模块包括坐标变换子模块、轴面流线计算子模块、轴面速度计算子模块、S1流面计算子模块、S2流面计算子模块、空间相对流线计算子模块、空间相对流线轴面速度计算子模块、空间相对水流角计算子模块、增厚系数试算子模块、虚拟厚度叶片排挤系数计算子模块和增厚系数插值子模块；所述前处理模块中的坐标变换子模块从所述数据采集模块获取所述全局参数和木模图参数；所述坐标变换子模块输出所述全局参数、木模图参数的轴面投影形式和木模图参数在S1流面坐标系下的表现形式，分别传输至轴面流线计算子模块、轴面速度计算子模块、S1流面计算子模块和S2流面计算子模块；所述轴面流线计算子模块、S1流面计算子模块和S2流面计算子模块输出各自模块的计算结果，都传至空间相对流线计算子模块；所述空间相对流线计算子模块的第一路输出连同所述轴面速度计算子模块的结果一起传输至空间相对流线轴面速度计算子模块，第二路输出传至空间相对水流角计算子模块，第三路输出连同增厚系数试算子模块的结果一起传输至虚拟厚度叶片排挤系数计算子模块；增厚系数插值子模块从数值计算模块获取初始输入值，对已求得的导叶开度-增厚系数进行插值，计算结果返回到虚拟厚度叶片排挤系数计算子模块；所述前处理模块用于计算得到求解目标参数所需要的准备参数，并将所述准备参数传输至所述显示模块的数据存储子模块；

所述数值计算模块包括数值计算子模块和计算结果验证子模块，其中，所述数值计算子模块的输入是来自空间相对流线计算子模块、空间相对流线轴面速度计算子模块、空间相对水流角计算子模块和虚拟厚度叶片排挤系数计算子模块的准备参数；所述数值计算子模块的输出为目标参数的计算值和增厚系数；所述目标参数的计算值输出至计算结果验证子模块；所述计算结果验证子模块通过CFD数值模拟求得所述目标参数的模拟值，并且求得所述目标参数的计算值和模拟值之间的相对误差；若所述相对误差超出判定标准值，则改变增厚系数试算子模块中增厚系数的初始赋值，返回到虚拟厚度叶片排挤系数计算子模块并再次计算；当相对误差小于判定标准值，即认为所述目标参数的计算值是最终解析解，且所述增厚系数满足误差精度要求；选取一定数量的离散工况，每隔一个工况执行上述过程，将每个工况的增厚系数和对应的导叶开度传递至所述前处理模块的增厚系数插值子模块，插值得到未执行上述过程工况对应的增厚系数，再依次执行上述子模块，得其他导叶开度下所述目标参数计算值和相对误差；将所述计算值、增厚系数、CFD模拟值和相对误差一路传输至所述显示模块的计算结果显示子模块进行显示和输出，另一路传输至所述显示模块的数据存储子模块进行存储；

所述显示模块包括计算结果显示子模块和数据存储子模块，所述计算结果显示子模块用于显示和输出所述目标参数的最终计算值、所述增厚系数、所述CFD模拟值和所述相对误差，所述数据存储子模块用于存储所述准备参数、所述目标参数计算值、所述增厚系数、所述CFD模拟值和所述相对误差。

7.根据权利要求6所述的一种基于虚拟厚度的反击式转轮全局特性数值计算系统，其特征在于：所述空间相对流线上各节点处相对水流角计算基于流体相对流动方向与所述空间相对流线切向方向一致的假设，计算结果为所述空间相对流线与水平线之间夹角的组合插值。

Images