CFD风速场重建详解Word格式文档下载.docx

CFD风速场重建详解Word格式文档下载.docx

《CFD风速场重建详解Word格式文档下载.docx》由会员分享，可在线阅读，更多相关《CFD风速场重建详解Word格式文档下载.docx（12页珍藏版）》请在冰豆网上搜索。

从质量、动量以及能量守恒的基本原理出发构建流体流动的基本方程。

设定流场边界条件。

用ANSYSICEM15.0预处理软件建立地风场几何模型,在计算区域内划分合理的网格。

采用基于有限体积法的ANSYSFlUENT15.0软件，设定计算区域边界条件，采用数值模拟计算方法，得出整个风场网格节点数据、风流动的速度云图、速度矢量图以及速度流线图等等，在CFD-post后处理软件中进行流场数据后处理。

利用matalab软件计算出网格节点风速的梯度，找出速度梯度的极值点，作为布置传感器的测点，风场重建主要基于ECT图像重构思想，找出测点速度与全场速度的关系，进而求出测量矩阵，达到重建风速场的目的。

2风场的CFD模拟计算

2.1模拟风场前预处理

2.1.1模型介绍

利用CFD软件模拟复杂地形时，首先要考虑的是计算域的尺度。

由于实际地形，连续群山之间的风场流动的影响是相互的，在起伏不定的山丘地形中，由于受到前方山坡风场流动的影响，后方山坡也不能仅仅看成孤立山来模拟。

同理，后方山坡的气流影响，也会波及到前方山坡风速和湍动能的分布，因而在设定模拟的山地目标模型时，就需要设置合适的计算区域，才能保证研究目标不受其它山包风速流动的影响，因此就需要在山坡前后延伸一定的尺度大小，因此计算域的选取是进行复杂地形风场模拟时首先考虑的问题。

查阅文献可知，山丘地形的轮廓通常符合以下模型，比如贝尔模型（Bellshape），余弦模型（Cosinesquared），高斯模型（Gaussian），正弦模型（Sinusoidal）等。

[1]

采用ANSYS15.0软件包里的ICEM15.0预处理软件构建合理的几何模型，在本次二维模拟中，两个小山包的轮廓线都符合正弦曲线的形式，第一个山包的正弦模型是:

（2-1）

第二个山包的正弦模型是：

（2-2）

其中，x表示山脊上各点方向的坐标值，m；

y表示山脊上各点的高度，m.

本文中，计算域最左下方的点设定为坐标原点（0,0），第一个山包的山底两个点的坐标分别为（25.28,0）、（31.58，0）；

第二个山包轮廓线与x轴的两个交点分别是（56.86,0）、（63.14,0）其中x是横坐标，计算区域两个山包的高度分别1m和2m,为了更好的模拟出在山包前后流体的流动，在山包前后分别将计算域水平延长了30m左右，两个山包山顶之间的距离为30m，总的计算域水平长度为100m，高度是30m，如图2-1所示：

图2-1计算区域几何形状

2.1.2网格划分

由于本论文考虑的是平面流动问题，所以在ICEM软件里需要生成面网格，对计算域进行网格划分，为了便于生成高质量的结构化网格，借助block,将计算域划分为八个块，采用的就是传统的映射网格划分技术，添加网格辅助线后创建的块如图2-2所示。

图2-2分块后的模型

在生成块之后，进行线的关联，然后调整块上点的位置，力求生成质量比较高的网格，为了得到山包附近更为细致的流动，需要在山包附近进行网格加密，本论文采用的是在块10沿x方向布置10个节点，而在上文我们知道第一个山包底面的长度是（m）,平均0.314m一个节点，而在远离山包的计算区域网格节点间隔就会大一些。

图2-3显示了网格生成结构。

图2-3生成的网格示意图

2.2计算域求解

2.2.1模型设置

本文采用ANSYS15.0软件包里的FLUENT软件，该软件采用有限体积法求解湍流Navier—Stokes方程，设定工作介质为空气，采用二维定常、不可压缩雷诺时均N-S方程描述计算域内部流场；

选用压力速度修正方法作为压力和速度的耦合方式，即SIMPLE算法，分析类型为稳态，计算收敛标准设置为1e-8。

由于现在只考虑二维平面，因此需要对模型选择ANSYSFLUENT二维求解器，导入FLUENT15.0里进行求解。

FLUENT15.0能够提供很多种湍流模型，在保证符合物理现象、节省时间，高精度的基础上，选择符合实际流体流动特点的湍流模型进行求解。

在本论文中选用带有壁面函数的标准湍流模型，一方程模型，标准的二方程模型分别来对模型进行求解。

确定最合适的求解湍流模型之后，再用此湍流模型对应求解结果作速度梯度运算，找出最适合布置传感器的测点。

入口风速设置为10m/s.

2.2.2计算结果对比

分别采用三种不同的湍流模型计算得到的速度分布，图2-4显示为标准二方程湍流模型的速度云图分布。

图2-4模型速度云图示意图

图2-5显示了采用标准二方程湍流模型的计算所得的计算域风场速度分布云图。

图2-5模型速度云图示意图

图2-6显示的了采用一方程湍流模型计算所得的风速在计算域的分布云图。

图2-6模型速度云图示意图

由以上不同湍流模型对应的速度矢量图和速度云图可以看出，改变了求解的湍流模型之后，流场速度分布云图，流场速度矢量图并没有明显的变化，由速度矢量图可以看到，在第二个山包后方都会形成回流。

2.2.3速度梯度的求解

在上一节中，我们求解得到了利用一方程模型，标准的二方程模型，标准的二方程模型模拟出来的风速场结果。

在速度分布云图中观察可知，它们并没有太大的区别，由于标准的k-e模型使用更为广泛，普适性更好，因而在本例中我们可以选用带有壁面函数的标准k-e湍流模型模拟风场。

节点数总共12144个，因此我们有12144个数据。

我们需要找到风场中最具代表性的点来作为传感器的布置位置。

如何找到这些特定点对接下来风场重建有着很重要的意义。

在ECT重建算法中我们有如下的公式：

（2-3）

在公式中，矩阵S是我们要求解的测量矩阵，一维列向量D中的元素由CFD模拟出来的网格节点数据组成，而传感器测点的速度值组成了一维向量T。

我们可以通过编制MATLAB语言计算得出网格节点的速度梯度，由于本文中在ICEM软件里设置的网格为结构化矩形网格，因此每个节点周围都包围着八个节点，但是他们不一定就是按照均匀矩形布置的，这是因为在山包附近的结构化网格如图2-7所示，而计算域中离山包较远的结构化网格由于不受山脊形状的影响，我们可以设置为如图2-8的均匀矩形网格。

图2-7山包周围的网格划分

图2-8远离山包区域的网格划分

由于速度梯度是沿着风向计算得出的，因而我们需要考虑每个节点速度的方向，由于附录一给出了每个节点的风速分别在x、y方向的矢量值以及每个节点的坐标，计算思路大致如下：

已知某个节点（x0,y0）的速度矢量是（u0,v0），可以求得每个节点速度方向与x轴的夹角，由于我们用FLUENT导出了每个节点的风速数据，我们可以通过比较风向与x轴的夹角与（450*n）角的关系，其中n=0，1，2，3，4，5，6，7，8。

风向角与（450*n）角最接近的对应角度，我们将其赋值给风向角，即为测点风向角，这样处理就方便我们找出与测点在同一个风向上的前一个节点，相当于将测点后移，这样的点找到之后，然后利用公式2-4求解梯度值。

（2-4）

其中是梯度值，是两个节点之间的距离。

要求解出各个网格节点的梯度值我们需要利用MATLAB编制程序来实现。

由于计算域入口边缘的网格节点上的速度梯度是求不出来的，而且这些节点对整个风场的意义不大，因而可以去除掉这些节点，这样我们就可以得到68*175个节点的速度梯度值。

同时我们也在MATLAB软件中得到了沿着各个网格节点的速度梯度分布曲线，如图2-9所示：

图2-9速度梯度分布

图2-9中，Z坐标代表着速度梯度值，在图中能够比较形象的看到那些速度梯度绝对值比较大的点分布，结合MATLAB计算出的68*175个梯度数据，从而确定它们的坐标值，即为第二种方案布置传感器的位置点。

从速度梯度分布曲线上可以看出，速度梯度比较大的地方主要在两个山包附近，根据MATLAB输出的风速梯度值，我们可以看出，在迎风坡前一段距离风速梯度值为负值，在接近山顶的位置处梯度值为正值且比较大，在背风坡一段区域风速梯度值为负值，结合图2-5显示的速度云图也可以看出背风坡一小部分的区域风速沿着风向变小了。

2.2.4均匀布置测点与按梯度极值点布置的比较

采用均匀布置测点，我们需要利用这些节点的数据来检测随着入口速度的变化，测点的速度变化情况。

具体做法就是，在同一个山地模型下，改变入口边界的初始速度为15m/s，利用相同的计算方法和湍流模型我们分别能得到风场网格节点速度数据和速度梯度数据。

我们认为风速变化的大小可以反映测点的灵敏度好坏，因而需要找到这些测点的速度平均变化值，可以求得每个测点速度变化的数值，再取其和求平均值，通过计算可得出测点均匀布置时，测点平均速度变化大小为5.065m/s。

同样的，求取每个测点在入口速度改变的情况下速度的变化值，再取其和求平均值，这一过程可以在通过在excle中计算实现，计算得出测点平均速度变化大小为11.081m/s。

通过以上两种方式，可以看出判断两种方式布置传感器的灵敏度，显而易见由于第二种布置测点的方式测点速度平均变化值比第一种方式更大，因而第二种测点布置方式即按照梯度分布布置更能反映风场的变化情况，这些点对入口风速的变化反映更为灵敏，总结出，通过按照速度梯度极值点来布置传感器的位置能够更好的反映整个风场的速度，这些测点的值即构成了公式中的中的元素。

2.3水平方向的二维风场计算

2.3.1几何建模

在前文中我们模拟了竖直方向的二维山脊，探索了均匀布置传感器与依据风速梯度布置传感器的优劣。

由于在竖直平面，我们无法改变入口风速的风向，为了研究速度梯度与入口风速风向的关系，在本节中，我们仍以前文的正弦山包模型为例，作出其在山高为0.5m的水平面投影，山包的投影分别是半径为1.047m和2.636m的圆。

x方向

图2-10水平方向二维计算域示意图

图2-10显示了水平方向二维山包计算域，其中计算域长度为100m，宽度为70m。

山后方计算域加长了50m左右，前方计算域加长20m左右，由于风从山前方吹过来，因此山前方计算域不需要加很长，两座山之间的距离仍然保持和前文一致。

2.3.2计算域求解

我们选择带有壁面函数的标准湍流模型，计算域左边界和下边界设置入口边界条件velocity-inlet，定常来流合速度设置为10m/s，在本例中依次设置入口风向角与x方向的夹角为0度，30度，60度；

计算域上边界和右边界均设置为出口边界条件outflow；

山包壁面边界条件设置为无滑移壁面wall。

设定工作介质为空气，采用二维不可压缩雷诺时均N-S方程描述计算域内部流场；

选用压力速度修正方法作为压力和速度的耦合方式，即SIMPLE算法，分析类型为瞬态，计算收敛标准设置为1e-8。

可在FLUENT中求解出的入口风速与x方向成0度，30度以及60度时的某同一时刻的计算域风速云图，每组分别输出了21422个网格节点的数据，这些数据可以用于接下来风速梯度的求解。

2.3.3不同入口风速方向的风速梯度对比

通过编制MATLAB语言，可以实现网格节点的风速梯度求解。

首先我们要将CFD-post里输出的节点数据进行规则化，即插值在规则的网格