CFD中文网

coolhhh

@BznW 压力脉动问题主要是因为入口质量通量不平衡导致的，NSRFG原论文中算出来的脉动压力系数也是非常大。简单的入口质量通量修正可以参考以下文章：
1.《Kim, Y., I.P. Castro and Z. Xie, Divergence-free turbulence inflow conditions for large-eddy simulations with incompressible flow solvers. Computers & Fluids, 2013. 84: p. 56-68.》
2.《Wang, Y. and X. Chen, Simulation of approaching boundary layer flow and wind loads on high-rise buildings by wall-modeled LES. Journal of Wind Engineering and Industrial Aerodynamics, 2020. 207: p. 104410.》

coolhhh

@李东岳那理解了，谢谢李老师

coolhhh

@李东岳李老师，看了boundaryFoam解析，也提到了充分发展湍流特性，但我还是没太明白充分发展湍流为何是这种特性，主要疑问如下：

对于充分发展的槽道湍流，瞬时的$\frac{\partial U(t)}{\partial x}$应该不等于0（因为流场为湍流，同个高度速度都有脉动值），这里的特征$\frac{\partial U}{\partial x}=0$中的速度是否指平均速度？

对于不可压缩流体，入口只施加平均风剖面速度（比如$U_x(z)$用对数率剖面，$U_y(z)=0$，$U_z(z)=0$）。之前算过这种工况的LES，结果是整个流场平均均按照自保持发展，基本没改变。对应的数学描述为：$\frac{\partial U_x(t)}{\partial x}=0$，$\frac{\partial p}{\partial x}=很小值$。另外再根据知乎帖子，这种平均风流动是否可称为充分发展剪切层流？

知乎：流体力学自由来流、自由剪切湍流、自由剪切流？

之前对湍流充分发展概念不理解，然后自己做了个猜想。把湍流发展类比于高中物理中的经典两个小球非弹性碰撞问题，湍流发展发生碰撞，碰撞就会产生力，力改变速度，接着发展趋向于使得两个碰撞物体平均速度（但瞬时速度不同）变为一致，最后脉动速度以相同平均速度向前发展形成充分发展湍流。
（1）最后脉动速度以相同平均速度向前发展形成充分发展湍流，这正对应于泰勒冻结假设：湍流在顺风向以平均速度移动，且湍流在运动过程中没有发生形态的显著变化。
（2）针对上面阐述，很容易产生疑问：充分发展时，湍流微元体在同个高度不同位置的速度不同，但为何都只是按照平均速度在输运？脉动速度对输运不起作用吗？针对这个问题尝试做的解释如下：
①x向的动量如下，大写字母表示合速度，小写表示脉动速度，脉动速度均值等于0

②将式(2)代入式(1)，并以移项得到式(3)：

③当湍流充分发展，微元体受到的合力应该等于0（不引起速度改变），即式(3)右边项等于0。根据泰勒冻结假设，湍流形态没显著变化，对应剪切应力不变，因此充分发展的湍流场：$\frac{\partial p}{\partial x}=常数$，并由于空气的粘性很小，对应剪切力很小，因此$\frac{\partial p}{\partial x}=很小值$。这对应了风工程中研究大气边界层湍流风场生成方法，为何要求空风洞中的压力接近于很小值才合理。

④因为充分发展流场按照泰勒冻结假设发展，这正对应于式(3)左边括号第一项等于0，即脉动速度物理量以平均速度输运：

⑤那么最后式(3)左边括号第2项也应等于0。该项等于0，说明了不同脉动速度的输运会相互抵消，这解释了为何充分发展时，不同位置的瞬时速度不同，但都以平均速度输运。

以上是困扰我很久的问题，针对自己疑问做的猜想，不知道合不合理，请各位大佬批评指正。总结一下猜想：三维充分发展湍流的特征是按照泰勒冻结假设方式发展，流场中压力梯度为很小值；瞬时脉动速度梯度不等于0：$\frac{\partial u}{\partial x}\neq 0, \frac{\partial u}{\partial y}\neq 0, \frac{\partial u}{\partial z}\neq 0 $

coolhhh

@李东岳或许是我这边fluent有某个设置错了，杨老师论文中顶面用的slip结果不会出现k衰减。我再试了fluent中顶面设置no-slip，结果是比较符合预期，接近顶面壁面平均速度减小，湍动能增大
Ux剖面.png TKE剖面.png

coolhhh

似乎找到原因了，是要把算例0文件夹中的U、k、epsilon中用codedFixedValue编写的边界name命名为不同名称，基于李老师算例文件修改后如下，用OpenFOAM-9直接Allrun
Yang2009OFCase.tar

    inlet
    {
        type            codedFixedValue;
        value           uniform 0.1; //default value
        name            linearTBC1; //name of new BC type

修改后的算例结果如下，近壁面有点偏差，但总体是保持住了

coolhhh

@李东岳李老师，我尝试按照自己的理解，用fluent尝试算了Yang2009算例，参照OpenFOAM的设置，结果能对上，但和原论文结果有点区别：

按照Yang2009设置，顶面设置为slip wall，结果基本能自保持，但接近顶部的k无法自保持，猜测是顶面边界条件设置问题。

fluent中slip wall设置如下，把壁面剪切应力设置为0。

把顶面设置为symmetry，结果就基本一样。实在不知道为何OpenFOAM算出来k衰减那么多

coolhhh

@coolhhh 在中性大气环境湍流动能的自保持 | 附有算例下载中说：

尝试设置为simpleGrading (1 10.9 1)，网格总数不变，对应的网格和Yang2009中网格接近，第一层网格高度0.01m，增长率是1.06。计算的Y+数值在80左右，但计算结果还是跟上面的计算结果基本一样，k偏差仍较大。

这里更正，计算的Y+数值大概在160左右

coolhhh

@李东岳李老师，算了您给的算例，计算结果跟您的结果一样，还没找到原因。但有以下问题：

blockMeshDict中把inlet和outlet设置成了wall，但都改为patch后，结果还是跟上面的计算结果基本一样，k偏差仍较大。

算例中blockMeshDict初始设置

convertToMeters 1;

vertices
(
    (0  0   0)
    (12 0   0)
    (12 1.8 0)
    (0  1.8 0)
    (0  0   1.8)
    (12 0   1.8)
    (12 1.8 1.8)
    (0  1.8 1.8)
);

blocks
(
    hex (0 1 2 3 4 5 6 7) (120 42 18) simpleGrading (1 4 1)
);

edges
(
);

boundary
(
    top
    {
        type patch;
        faces
        (
            (3 7 6 2)
        );
    }
    inlet
    {
        type wall;
        faces
        (
            (0 4 7 3)
        );
    }
    outlet
    {
        type wall;
        faces
        (
            (2 6 5 1)
        );
    }
    bottom
    {
        type wall;
        faces
        (
            (1 5 4 0)
        );
    }
    sides
    {
        type wall;
        faces
        (
            (0 3 2 1)
            (4 5 6 7)
        );
    }
);

Yang2009中网格，第一层网格高度0.01m，增长率是1.06。simpleGrading (1 4 1)这个设置反算对应第一层网格高度为0.0197，网格增长率为1.034，计算Y+结果大约在300：

# y+ ()       
# Time        	patch         	min           	max           	average       
0             	bottom	3.602087e+02	3.602087e+02	3.602087e+02
0             	sides	0.000000e+00	0.000000e+00	0.000000e+00
100           	bottom	3.026569e+02	4.241725e+02	3.196787e+02
100           	sides	2.660222e-03	1.234693e+02	2.231677e+01
200           	bottom	2.978864e+02	4.241725e+02	3.151323e+02
200           	sides	2.712471e-03	1.696112e+01	7.879751e-01
300           	bottom	2.978864e+02	4.241725e+02	3.151323e+02
300           	sides	3.569842e-03	1.696111e+01	5.413341e-01
400           	bottom	2.978864e+02	4.241725e+02	3.151323e+02
400           	sides	1.932649e-03	2.113821e+01	6.643026e-01
500           	bottom	2.978864e+02	4.241725e+02	3.151323e+02
500           	sides	1.392717e-03	1.696111e+01	5.073163e-01

尝试设置为simpleGrading (1 10.9 1)，网格总数不变，对应的网格和Yang2009中网格接近，第一层网格高度0.01m，增长率是1.06。计算的Y+数值在80左右，但计算结果还是跟上面的计算结果基本一样，k偏差仍较大。


# Time        	patch         	min           	max           	average       
0             	bottom	1.868099e+02	1.868099e+02	1.868099e+02
0             	sides	0.000000e+00	0.000000e+00	0.000000e+00
100           	bottom	1.598811e+02	2.174476e+02	1.701706e+02
100           	sides	2.192814e-02	1.108888e+02	2.684494e+01
200           	bottom	1.539062e+02	2.174475e+02	1.625827e+02
200           	sides	8.207370e-04	2.373603e+01	1.893065e+00
300           	bottom	1.537618e+02	2.174475e+02	1.625742e+02
300           	sides	4.922483e-03	8.914920e+00	7.251943e-01
400           	bottom	1.537618e+02	2.174475e+02	1.625742e+02
400           	sides	1.176434e-03	3.476801e+01	9.170413e-01
500           	bottom	1.537618e+02	2.174475e+02	1.625742e+02
500           	sides	3.049010e-03	8.458892e+00	7.410303e-01

Yang2009中，提到了SKE1和SKE2两种模型，并模拟对比。

尝试根据式(25)和(26)把这两个k曲线画出对比，这两条曲线还是有挺大区别，但在原论文中两个工况的inlet结果看起来像是一样的

coolhhh

@李东岳收到，谢谢老师

coolhhh

@李东岳李老师，这是我的邮箱 coolclw67@163.com。之前没研究RANS，但挺感兴趣的，最近我也研究下，谢谢李老师

coolhhh

@coolhhh

coolhhh 发布的最新帖子