天然气管道运行模拟及仿真技术研究.doc

天然气管道运行模拟及仿真技术研究
1011202045 蔡永军 科学计算选讲结课论文
为了预测天然气管道运行状态，制定合理的管输计划，更好的配置设备开机，天然气管道输送过程中需要进行工况模拟及仿真。实际工作中需要建立压缩机、阀门等设备的模型，确定管段的控制方程、气体的状态方程，针对给出的初始条件和边界条件，筛选确定天然气管网数学模型的离散方法与非线性方程组的求解算法寻找合理的非线性方程的求解算法，得到合理的数值解。

1天然气管道仿真数学模型
1.1管段的控制方程
对于管道中的任意管段，经过适当的简化可以用下列公式来描述：
连续性方程：
（1）
运动方程：
（2）
能量方程：
（3）
式中：A——管道的横截面积，m2；
ρ——流体密度，kg/m3；
t——时间，s；
x——坐标，m；
u——速度，m/s；
P——压力，Pa；
θ——管道倾角，rad；
λ——水力摩阻系数；
D——管道内径，m；
T——流体温度，k；
k1——流体至管壁的换热系数；
h——比焓；
Tw——管壁的温度，k。
1.2 阀门控制方程
阀门控制方程如下：
（4）
式中： Mup——阀门入口质量流量，kg/s；
Mdw——阀门入口质量流量，kg/s
Cg——阀门系数；
Pup——阀的入口压力，Pa；
Pdw——阀的出口压力，Pa。
1.3压缩机控制方程
简化后的压缩机控制方程如下
（5）
式中：——压缩机压比；
m——多变压缩指数；
n——压缩机的实际转速，rpm；
n0——压缩机的额定转速，rpm；
a1, b1, c1——系数；
Q——给定状态下的体积流量，m3/s；
1.4 理想调节阀阀控制方程
理想调节阀控制方程如下：
（6）
2气体的状态方程
采用BWRS气体状态方程，如下：
（7）
式中：P——系统压力，KPa；
T——系统温度，K；
ρ——混合气体密度，Kmol/m3；
R——气体常数，8.3143KJ/( Kmol.K)。
A0、 B0、C0、D0、E0、a、b、c、d、α、γ为方程的是一个参数，根据（8）确定。
（8）
式中：xi、xj——混合气体中i和j组分的摩尔分数；
kij——为i、j组分间的交互作用系数。
3气体的焓方程
气体的焓方程如下：
（9）
4 管道周边的热力模型
管道的有效土壤厚度采用等效圆筒法，传热半径由下式计算：
（10）
式中：R2-R1——土壤厚度，m；
R1——从管道中心至土壤层的半径，m；
H——至管道中心的实际埋深，m；
D——管道直径，m。
管道和周围环境的瞬态热力模型计算式如下：
（11）
式中：k——周围环境导热系数；
r——传热半径；
Tr——r处的气体温度；
Cp——气体定压比热；
Tt——t时刻的气体温度。
单位管长热流量由下式表示。通过该公式计算管壁在任意节点的温度。
（12）
式中：k2——管壁至土壤换热系数；
K1——流体至管壁换热系数；
Tw——管壁温度；
T0——R2处的温度；
T——气体温度。
5 水力摩阻系数计算式
管段控制方程涉及的水力摩阻系数λ采用F.Colebrook－White公式计算，该公式表达如下：
（13）
式中，e/D——管道粗糙度和内径的无因次比；
Re——雷诺数。
6控制方程的离散化
由管道控制方程与气体状态方程组成的非线性偏微分方程组，一般不能得出管流气体基本变量的解析解，因此有必要应用计算数学的方法求解偏微分方程组的数值解。本专题中选用中心隐式差分法对控制方程进行离散化。确定采用的基本变量为气体的密度（ρ）、速度（u）和温度（T）。
6.1离散形式
引进变量φ，φ代表三个流动基本中的任意一个。在时间步长为Δx , 空间步长为Δt 的情况下，以空间i和时间网格点t采用中心隐式差分格式，则有以下离散形式：
对于基本流动变量：
（14）
基本流动变量对时间的一阶偏导数：
（15）
基本流动变量对空间的一阶偏导数：
（16）
基本流动变量对时间的二阶偏导数：
（17）
基本流动变量对空间的二阶偏导数：
（18）
基本流动变量对空间及时间的二阶偏导数：
（19）
6.2 离散后的控制方程
离散后的控制方程如下：
离散后的连续性方程：
（20）
离散后的运动方程：
（21）
离散后的能量方程：
（22）
6.3 初始条件与边界条件
初始条件指系统开始运行时的初始压力、流量或温度的分布状态。边界条件指某一管段起始节点和终止节点上的约束条件。主要包括：
(1) 管段端点上的输油泵、压缩机或阀门等的出入口压力、流量、温度、转速、压比或开度设定值；
(2) 气源对应节点的压力、流量或温度设定值；
(3) 分输点对应节点的压力、流量或温度设定值；
(4) 节点处压力、流量或温度的一致性；
(5) 节点处压力、流量或温度的范围控制值；
(6) 管道物理元件周围的温度场状况。
7非线性方程组的求解算法
离散后的控制方程配合边界条件和初始条件才能封闭，封闭后形成了非线性方程组，对于该非线性方程组选取牛顿迭代法进行求解。
若采用的矩阵形式（其中C(x)为非线性方程组的系数矩阵），则为需要求解的向量，为等式右边的向量。
（1） 牛顿拉普森迭代法
设迭代函数列
迭代变量
迭代增量
迭代函数
牛顿拉普森迭代公式如下：
（5.7-1）
对于迭代函数F， 将求解非线性方程组问题转化成为寻根问题，也即要求下式成立：
（5.7-2）
对任意点x0和它的相邻点/邻域（x0+△x）， 通过泰勒展开式我们有：
（5.7-3）
若采用矩阵形式，则有：
（5.7-4）
其中 J 为n×n的雅可比矩阵且。
如果略去其中的高次项，并要求。我们得出：
（5.7-5）
至此，可以按照牛顿拉普森迭代法的求解步骤进行计算。
（2） 牛顿＋最速搜索迭代法
引入目标函数 ，通过简单的数学运算，可以得到牛顿迭代沿着此目标函数的梯度方向，始终可以发现一个数值α能使得目标函数的值下降，即：。所以本专题研究采用如下迭代公式：
（5.7-6）
该方法较牛顿拉普森方法具有收敛速度快，且全局收敛的特点。
8仿真运行
8.1仿真运算的基本过程
仿真运行的基本过程如图1仿真运行所示。运行系统状态数据和不同视角构成单文档－多视的关系。
图中的兰色带箭头线条表示通过不同的视角和核心功能层接口，可以监视或编辑仿真系统的状态数据，并下达计算命令；粉色线条表示运算中仿真模拟器和数据模块进行数据交互；绿色线条表示初始化过程加载数据。

图1仿真运行过程
仿真运行基本过程如下：
(1) 通过各个接口或视角，下达加载系统命令；
(2) 调度模块命令编译模块加载编连文件；
(3) 如果需要，通过各个接口或视角，对系统参数进行进一步初始化；
(4) 通过各个接口或视角，下达一轮计算命令；
(5) 仿真模拟器开始一轮计算，并输出结果到数据模块；调度模块通知相应接口和视角计算结束；
(6) 相应接口或视角获取关心的数据；
(7) 如果需要进行新一轮计算，相应接口和视角可以对部分参数进行重新设定，并下达新一轮计算命令，系统将回到第5步。
8.2仿真运算的实现结构
仿真运算由仿真模拟器作为核心模块来实现，两个直接的辅助模块是数据模块和编译模块。
这3个模块相互协同进行仿真运算，基本过程如图2所示。

图2仿真模拟器结构
(1) 编译模块根据模型文件和元件库，建立编连文件，编连文件中包含仿真系统各个元件的信息；同时，根据通用、基本的控制方程组、元件参数关系计算函数库、外部边界方程组，编译模块将建立此系统的全部方程组，一并加入到编连文件中；
(2) 在运算前的初始化过程中，编连文件中的参数数据将被加载到数据模块，方程组将被加载到仿真模拟器相应的方程组列表中，包括：控制方程组列表、外部边界方程组列表、元件参数关系计算函数库列表；
(3) 每轮运算开始前，仿真模拟器的输入参数预处理模块需要对输入参数进行预处理，例如对部分参数进行离散化或拟合；
(4) 开始运算后，仿真模拟器的数值计算模块根据处理好的输入参数和相关的方程组进行数值计算；
(5) 数值计算模块计算出的是需要联立求解的基本输出参数；
(6) 根据输入参数、基本输出参数、元件参数关系计算函数库，仿真模拟器同步对其它输出参数进行数值化计算；
(7) 其它输出参数被输出到数据模块，此后可以通知调度模块前来获取相关数据并进行下一轮计算的参数输入了。
9结论
通过建立天然气管道数学模型及求解算饭，在输入管道的基本参数后，可以根据输入的初始运行状态预测下一时刻的运行状态，从而为排定管输计划、优化运行工况提供决策依据。