一种管道瞬变流的计算方法

1.本发明涉及能源动力工程仿真技术领域，具体而言，涉及一种管道瞬变流的计算方法。


背景技术：

2.管道瞬变流也常称为水击，水击过程时常伴随着巨大的噪声和剧烈振动，严重时会造成管道疲劳断裂以及设备故障。常用的管道瞬变流的计算方法有基于欧拉格式的特征线法、有限体积法、有限差分法等，这类计算方法是目前常用的瞬变流计算方法，但是存在一个普遍的问题，就是复杂管道瞬变流计算效率低下，难以达到实时计算。另一类是基于拉格朗日格式的波特性法、拉格朗日法、波追踪法等，这类方法最大的特点就是追踪压力波的传播过程，并在该波经过的地方更新流场信息，因此具有计算效率高的优点。
3.理论上看，基于拉格朗日格式的数值方法具有计算效率高的优势，但是受制于计算方法的实施方案不成熟，目前这类方法的计算效率相对欧拉格式的数值方法的效率优势并不明显，存在压力波数量指数增长、压力波消失阈值的选取等问题。这些问题会导致计算精度和计算效率低，在牺牲一定计算精度的前提下最多只能提升60％的计算效率。采用准稳态摩阻或者恒定流摩擦进行瞬变流计算也会导致计算精度低。
4.目前的瞬变流计算技术方案主要是基于特征线法，导致为提高计算效率必须舍弃部分计算精度，比如在管道节点划分时采用较大的波速调整系数和短管划分比例系数，虽然能够减少管道节点划分总量，但会引入不可忽视的计算误差。
5.因此有必要提供一种新的管道瞬变流的计算方法，以克服现有的管道瞬变流的计算方法所存在的计算效率低或者计算精度低的问题。


技术实现要素：

6.本发明提供一种管道瞬变流的计算方法，用以克服现有技术中存在的至少一个技术问题。
7.根据本发明实施例技术方案，提供一种管道瞬变流的计算方法，包括：
8.步骤s1、设计基于非恒定流阻力模型的高效波追踪方程组，所述高效波追踪方程组的表达式具体包括：
[0009][0010][0011]
其中，符号fa表示所述管道的端点a处的下游波，fb表示所述管道的端点b处的下游波，符号l
ab
表示管道的端点a和所述管道的端点b间的距离；符号g表示重力加速度；符号j表示所述管道的当前摩擦损失，符号j0表示所述管道的初始摩擦损失；下角标a、b表示单根管道的上下游端点，(j-j0)a表示所述端点a到所述端点b的压力波损失，(j-j0)b表示端点b到
端点a的压力波损失，式(1)和式(2)中的各项可由下式计算得到：
[0012][0013][0014][0015]
其中，符号f表示初始流速正方向的下游波；符号f表示初始流速负方向的上游波；符号h表示静压头，符号q表示体积流量；符号a表示波速；符号d表示所述管道的直径；符号ac表示所述管道的横截面积；符号λ表示所述管道的沿程阻力系数，符号c
*
表示取值与流态有关的brunone阻力系数；符号q0表示管道稳态时的流速；符号qa表示所述端点a处的流速，符号qb所述端点b处的流速，符号λa表示所述端点a处的沿程阻力系数，符号λb表示所述端点b处的沿程阻力系数，q
a,0
所述端点a处的初始稳态流速，q
b,0
所述端点b处的初始稳态流速，符号λ
a,0
表示所述端点a处的初始稳态沿程阻力系数，符号λ
b,0
表示所述端点b处的初始稳态沿程阻力系数；
[0016]
步骤s2、确定短管划分标准，进行管道节点划分以及时间步长的确定；
[0017]
步骤s3：给定初始条件和边界条件，进行稳态计算；
[0018]
步骤s4：具体进行高效波追踪法的计算，具体包括：
[0019]
步骤s41、每一根管道都建立一个基于空间域的压力波库，压力波库由两个分别存储管道上游波f和下游波f的数组组成，数组大小为n+1；下游波数组中的元素0,1
…
n与管道空间节点位置一一对应，上游波数组元素与管道节点位置的对应关系相反，数组元素存储管道相应位置处的压力波幅值；压力波n表示传入管道端点的压力波，压力波0表示管道端点传出的压力波；
[0020]
步骤s42、在每一个时间步长下，调用传播到边界节点处的压力波n，进行边界节点方程的求解，由边界处产生的压力波与入射波的方向是相反的，上游边界节点a处进入的压力波为上游压力波库f中的压力波n，产生的压力波为下游压力波库中的压力波0，新产生的压力波0在插入到数组前先进行从节点a到节点b的压力波衰减计算，新产生的压力波记为f0，f0，由式(6)和式(7)来进行衰减计算：
[0021][0022][0023]
对于新产生的压力波0直接从数组元素0处插入数组中，使数组中的元素朝着压力波传播的方向前进一个位置；数组元素n被调用后直接删除；
[0024]
步骤s43、在计算过程中追踪压力波f和f，通过下式获取边界节点的压头和体积流量，更新管道端点的瞬时水力状态，在下一个时间步长下返回步骤s2，直至完成瞬态计算；
[0025][0026]
优选的，所述管道节点划分具体包括：
[0027]
步骤s21、将管道系统中的管道按照长到短排序；
[0028]
步骤s22、选择最长管道l
max
，根据公式l＜δl
max
确定短管长度l的划分标准，符号δ表示用以确定短管的长度标准的比例系数；步骤s23、在短管划分结束后，将管道系统中的短管剔除，在剩余管道中选取最短管道l
min
；
[0029]
步骤s24、令最短管道分段数n
min
＝1，根据最短管的波速计算时间步长
[0030]
步骤s25、采用步骤s24中的时间步长以及剩余各管的波速进行节点划分，当计算的分段数不是整数时，通过调整波速来获得整数分段数，记录每根管道的波速调整系数
[0031]
步骤s26、判断所有管道的波速调整系数是否小于最大波速调整系数若满足，则管道划分结束；否则，将最短管道分段数加1后，继续步骤s24，直至满足条件。
[0032]
优选的，所述初始条件包括管道各节点处的压力、流量、波速以及横截面积。
[0033]
优选的，所述边界条件包括但不限于阀门、水库。
[0034]
本说明书一个实施例至少能够达到以下有益效果：
[0035]
1.在进行管道节点划分中，短管划分比例系数和波速调整系数会影响计算精度和计算效率，两个系数越小，计算精度越高但计算效率将会降低。本发明中采用高效波追踪法能够解决该问题，因为高效波追踪法的计算时间对管道分段数不敏感，所以可以尽可能的满足精度需求而选取小的短管划分比例系数和波速调整系数。
[0036]
2.高效波追踪法的计算实施方案可以模拟随着时间推进的压力波在管道中的传播过程，还可以在每一个时刻进行管道的瞬时水力状态的更新，结合了拉格朗日格式和欧拉格式的优点，极大地提高了高效波追踪法的可操作性和计算高效性。
[0037]
3.本发明提出的基于高效波追踪法的管道瞬变流计算技术方案有效地提高了波追踪法在管道瞬变流计算上的效率和精度。根据本发明计算的案例中，具有与特征线法同等精度，计算效率较特征线法提升达90％以上。
附图说明
[0038]
为了更清楚地说明本发明实施例或现有技术中的技术方案，下面将对实施例或现有技术描述中所需要使用的附图作简单地介绍，显而易见地，下面描述中的附图仅仅是本发明的一些实施例，对于本领域普通技术人员来讲，在不付出创造性劳动的前提下，还可以根据这些附图获得其他的附图。
[0039]
图1是本发明提供的一种管道瞬变流的计算技术方法的流程图；
[0040]
图2是本发明技术方案提出的基于高效波追踪法的计算实施方案图；
[0041]
图3是本发明技术方案提出的基于高效波追踪法的管道瞬变流计算技术方案的计算结果与实验值的对比图；
[0042]
图4是本发明技术方案提出的基于高效波追踪法的管道瞬变流计算技术方案的计算效率与特征线法的对比图。
具体实施方式
[0043]
为使本说明书一个或多个实施例的目的、技术方案和优点更加清楚，下面将结合本说明书具体实施例及相应的附图对本说明书一个或多个实施例的技术方案进行清楚、完整地描述。显然，所描述的实施例仅是本说明书的一部分实施例，而不是全部的实施例。基于本说明书中的实施例，本领域普通技术人员在没有做出创造性劳动前提下所获得的所有其他实施例，都属于本说明书一个或多个实施例保护的范围。
[0044]
应当理解，尽管在本技术文件中可能采用术语第一、第二、第三等来描述各种信息，但这些信息不应限于这些术语。这些术语仅用来将同一类型的信息彼此区分开。
[0045]
如图1所示，本发明提供的基于高效波追踪法的管道瞬变流计算技术方案，包括如下步骤：
[0046]
步骤s1：采用基于非恒定流阻力模型的高效波追踪方程组；
[0047]
基于非恒定流阻力模型的高效波追踪方程组表达式为：
[0048][0049][0050]
其中
[0051][0052][0053][0054]
式中，符号fa表示所述管道的端点a处的下游波，fb表示所述管道的端点b处的下游波，符号l
ab
表示管道的端点a和所述管道的端点b间的距离；符号g表示重力加速度；符号j表示所述管道的当前摩擦损失，符号j0表示所述管道的初始摩擦损失；下角标a、b表示单根管道的上下游端点，(j-j0)a表示所述端点a到所述端点b的压力波损失，(j-j0)b表示端点b到端点a的压力波损失。符号fa表示所述管道的端点a处的下游波，fb表示所述管道的端点b处的下游波，符号l
ab
表示管道的端点a和所述管道的端点b间的距离；符号g表示重力加速度；符号j表示所述管道的当前摩擦损失，符号j0表示所述管道的初始摩擦损失；下角标a、b表示单根管道的上下游端点，(j-j0)a表示所述端点a到所述端点b的压力波损失，(j-j0)b表示端点b到端点a的压力波损失。
[0055]
步骤s2：确定短管划分标准，进行管道节点划分以及时间步长的确定。
[0056]
其中，节点划分流程如下：
[0057]
步骤s21、将管道系统中的管道按照长到短排序；
[0058]
步骤s22、选择最长管道l
max
，根据公式l＜δl
max
确定短管长度l的划分标准，符号δ表示用以确定短管的长度标准的比例系数；步骤s23、在短管划分结束后，将管道系统中的短管剔除，在剩余管道中选取最短管道l
min
；
[0059]
步骤s24、令最短管道分段数n
min
＝1，根据最短管的波速计算时间步长
[0060]
步骤s25、采用步骤s24中的时间步长以及剩余各管的波速进行节点划分，当计算的分段数不是整数时，通过调整波速来获得整数分段数，记录每根管道的波速调整系数
[0061]
步骤s26、判断所有管道的波速调整系数是否小于最大波速调整系数若满足，则管道划分结束；否则，将最短管道分段数加1后，继续步骤s24，直至满足条件。
[0062]
在进行管道节点划分中，短管划分比例系数和波速调整系数会影响计算精度和计算效率，两个系数越小，计算精度越高但计算效率将会降低。本发明中采用高效波追踪法能够解决该问题，因为高效波追踪法的计算时间对分段数不敏感，所以可以尽可能的满足精度需求而选取小的短管划分比例系数和波速调整系数。经计算，一般选取δ＝0.01，即可满足精度需求。
[0063]
步骤s3：给定初始条件和边界条件，进行稳态计算
[0064]
步骤s1中基于非恒定流阻力模型的高效波追踪方程组都与管道初始参数有关，因此在进行瞬态计算前需要根据初始条件和边界条件进行稳态计算。初始条件包含管道各节点处的压力、流量、波速以及横截面积；边界条件包含几种常用边界条件，如阀门，水库等。
[0065]
步骤s4：高效波追踪法的计算实施方案
[0066]
(1)每一根管道都建立一个基于空间域的压力波库，压力波库由两个分别存储管道上游波f和下游波f的数组组成，数组大小为n+1。下游波数组中的元素0,1
…
n与管道空间节点位置一一对应，上游波数组元素与管道节点位置的对应关系刚好相反，数组元素存储管道相应位置处的压力波幅值。压力波n表示传入管道端点的压力波，压力波0表示管道端点传出的压力波；
[0067]
(2)在每一个时间步长下，调用传播到边界节点处的压力波n，进行边界节点方程的求解。由边界处产生的压力波与入射波的方向是相反的，比如上游边界节点a处进入的压力波为上游压力波库f中的压力波n，产生的压力波为下游压力波库中的压力波0。新产生的压力波0在插入到数组前需要先进行从节点a到节点b的压力波衰减计算。新产生的压力波记为，f0，f0，那么衰减计算可由式(6)和式(7)来表示：
[0068][0069][0070]
对于新产生的压力波0直接从数组元素0处插入数组中，使数组中的元素朝着压力波传播的方向前进一个位置。数组元素n被调用后直接删除。这一过程正好描述了压力波在边界处的消失和产生以及在管道两端之间的传播与衰减过程；
[0071]
(3)在计算过程中只需要追踪压力波f和f，即，通过下式获取边界节点的压头和体积流量，更新管道端点的瞬时水力状态。在下一个时间步长下返回步骤(2)，直至完成瞬态计算。
[0072][0073]
上述技术方案可以模拟随着时间推进的压力波在管道中的传播过程，还可以在每一个时刻进行管道的瞬时水力状态的更新，结合了拉格朗日格式和欧拉格式的优点，极大地提高了高效波追踪法的可操作性和计算高效性。
[0074]
图3和图4分别是本发明的计算精度和计算效率的验证结果：
[0075]
实验数据来源于soares实验，该实验装置采用的是水库-管道-阀门系统。水库液位46m，管道长度15.22m，管道内径0.02m，管道壁厚0.001m，体积流量1.33e-4m3/s，关阀时间0.018s。由图3可以看到，计算的最大压头较实验值偏大2％左右，具有保守性。同时计算的压力波衰减曲线和实验值高度重合，表明本发明的计算精度较高。由图4中的左图可知，随着分段数的增加，特征线计算时间始终大于高效波追踪法的计算时间，而高效波追踪法的计算时间对分段数不敏感。从图4中的右图，计算效率提升率随着分段数增加而提高，在分段数小于10时，计算效率提升率变化明显，最低为25％；当分段数大于20时，高效波追踪法的计算时间效率提升达90％以上。
[0076]
本领域普通技术人员可以理解：附图只是一个实施例的示意图，附图中的模块或流程并不一定是实施本发明所必须的。
[0077]
最后应说明的是：以上实施例仅用以说明本发明的技术方案，而非对其限制；尽管参照前述实施例对本发明进行了详细的说明，本领域的普通技术人员应当理解：其依然可以对前述实施例所记载的技术方案进行修改，或者对其中部分技术特征进行等同替换；而这些修改或者替换，并不使相应技术方案的本质脱离本发明实施例技术方案的精神和范围。