钝体建筑扰流的大涡模拟
-
@hitsc30 在 钝体建筑扰流的大涡模拟 中说：
当前是按照tamura教授的设置帮边界条件、离散格式、求解器、亚格子模型全换了。就是采用的网格和他们的有所区别（他们采用的snappyhexmesh），准备搞个和他们一样的网格试一试。实在不行就只能放弃了
方便把建筑物的详细参数，以及对标的那篇实验分享一下吗？我跑一下看看。压力系数就是一个静压，如果这个都对不上，高阶统计量就能不可能了，那LES的作用就被完完全全浪费。 就算一个建筑物的话，这个已经属于10多年的范畴了，2006年就能跑DNS了。
https://www.sciencedirect.com/science/article/pii/S0142727X06000476 -
@hitsc30 在 钝体建筑扰流的大涡模拟 中说：
@gtian 您好，我研究的对象是基于东京工艺大学TPU数据库的1：1：4模型的模拟（在和风洞同样缩尺1：400后，模型的尺寸为0.1m0.1m0.4m。α为0.25，计算域按照cost732推荐为长8m宽4m高度为1.6m）。TPU数据库参见下述链接：
http://www.wind.arch.t-kougei.ac.jp/system/eng/contents/code/w_it
对标的就是该风压数据库提供数据处理得到的平均风压系数和脉动风压系数。当前我的平均风压系数一直偏小了0.1；迎风向缩尺后尺寸为0.1m，顺风向缩尺后尺寸为0.1m,竖向的尺寸为0.4m。指数选取为0.25
-
@gtian 有一篇类似的文章，该文章研究的模型是1：1：5模型。我研究的是1：1：4模型，在TPU数据库上实际有四百个个测点，每个面各有一百个测点。
下述链接是1：1：5模型的文件，链接如下，国内有人做过1：1：4模型，基于fluent平台进行的模拟，结果还未见刊。
https://www.sciencedirect.com/science/article/pii/S0167610520303202
雷诺数，在建筑高度处大约30w。计算高度处模型高度0.4m(对应实际高度160m)，风速11m/s。风洞数据库参见链接：
http://www.wind.arch.t-kougei.ac.jp/system/eng/contents/code/w_it
祝大家五四青年节快乐！ -
更新一下，跑了几天。
- 两个时间尺度，
timestep \delta T ， 控制Co < 2.
总运行时长： 第一阶段250个大涡turnover time ： 250 * (0.4（H） / 0.5 (friction velocity); 第二阶段100个大涡turnover time .
2.
4. Mesh stats points: 3520103 faces: 10396800 internal faces: 10234752 cells: 3438592 faces per cell: 6 boundary patches: 23 point zones: 0 face zones: 0 cell zones: 0 Overall number of cells of each type: hexahedra: 3438592 prisms: 0 wedges: 0 pyramids: 0 tet wedges: 0 tetrahedra: 0 polyhedra: 0 Checking topology... Boundary definition OK. Cell to face addressing OK. Point usage OK. Upper triangular ordering OK. Face vertices OK. Number of regions: 1 (OK). Checking patch topology for multiply connected surfaces... Patch Faces Points Surface topology bottomWall 49088 49552 ok (non-closed singly connected) topWall 49152 49601 ok (non-closed singly connected) sides1_half0 12288 12593 ok (non-closed singly connected) sides1_half1 12288 12593 ok (non-closed singly connected) sides2_half0 3072 3341 ok (non-closed singly connected) sides2_half1 3072 3341 ok (non-closed singly connected) sides3_half0 1536 1799 ok (non-closed singly connected) sides3_half1 1536 1799 ok (non-closed singly connected) sides4_half0 768 1028 ok (non-closed singly connected) sides4_half1 768 1028 ok (non-closed singly connected) sides5_half0 256 514 ok (non-closed singly connected) sides5_half1 256 514 ok (non-closed singly connected) out1_half0 9216 9457 ok (non-closed singly connected) in1_half1 9216 9457 ok (non-closed singly connected) out2_half0 2304 2509 ok (non-closed singly connected) in2_half1 2304 2509 ok (non-closed singly connected) out3_half0 1152 1351 ok (non-closed singly connected) in3_half1 1152 1351 ok (non-closed singly connected) out4_half0 576 772 ok (non-closed singly connected) in4_half1 576 772 ok (non-closed singly connected) out5_half0 192 386 ok (non-closed singly connected) in5_half1 192 386 ok (non-closed singly connected) building 1088 1105 ok (non-closed singly connected) Checking geometry... Overall domain bounding box (0 0 0) (3.2 1.6 2.4) Mesh (non-empty, non-wedge) directions (1 1 1) Mesh (non-empty) directions (1 1 1) Boundary openness (8.71455e-17 -4.51544e-16 -4.1172e-16) OK. Max cell openness = 1.73472e-16 OK. Max aspect ratio = 12.6864 OK. Minimum face area = 0.00015625. Maximum face area = 0.00198225. Face area magnitudes OK. Min volume = 1.95312e-06. Max volume = 2.47781e-05. Total volume = 12.284. Cell volumes OK. Mesh non-orthogonality Max: 0 average: 0 Non-orthogonality check OK. Face pyramids OK. Max skewness = 3.18096e-13 OK. Coupled point location match (average 0) OK. Mesh OK. End - 两个时间尺度，
-
/*--------------------------------*- C++ -*----------------------------------*\ | ========= | | | \\ / F ield | OpenFOAM: The Open Source CFD Toolbox | | \\ / O peration | Version: 2.4.0 | | \\ / A nd | Web: www.OpenFOAM.org | | \\/ M anipulation | | \*---------------------------------------------------------------------------*/ FoamFile { version 2.0; format ascii; class volScalarField; location "0"; object p; } // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // dimensions [ 0 2 -2 0 0 0 0 ]; internalField uniform 0; boundaryField { bottomWall { type zeroGradient; } topWall { type symmetryPlane; } sides1_half0 { type symmetryPlane; } sides2_half0 { type symmetryPlane; } sides3_half0 { type symmetryPlane; } sides4_half0 { type symmetryPlane; } sides5_half0 { type symmetryPlane; } out1_half0 { type fixedValue; value uniform 0; } out2_half0 { type fixedValue; value uniform 0; } out3_half0 { type fixedValue; value uniform 0; } out4_half0 { type fixedValue; value uniform 0; } out5_half0 { type fixedValue; value uniform 0; } sides5_half1 { type symmetryPlane; } sides4_half1 { type symmetryPlane; } sides3_half1 { type symmetryPlane; } sides2_half1 { type symmetryPlane; } sides1_half1 { type symmetryPlane; } in1_half1 { type zeroGradient; } in2_half1 { type zeroGradient; } in3_half1 { type zeroGradient; } in4_half1 { type zeroGradient; } in5_half1 { type zeroGradient; } building { type zeroGradient; } } // ************************************************************************* ///*--------------------------------*- C++ -*----------------------------------*\ | ========= | | | \\ / F ield | OpenFOAM: The Open Source CFD Toolbox | | \\ / O peration | Version: 2.4.0 | | \\ / A nd | Web: www.OpenFOAM.org | | \\/ M anipulation | | \*---------------------------------------------------------------------------*/ FoamFile { version 2.0; format ascii; class volVectorField; location "0"; object U; } // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // dimensions [ 0 1 -1 0 0 0 0 ]; internalField uniform ( 0.01335 0 0 ); boundaryField { bottomWall { type fixedValue; value uniform ( 0 0 0 ); } topWall { type symmetryPlane; } sides1_half0 { type symmetryPlane; } sides2_half0 { type symmetryPlane; } sides3_half0 { type symmetryPlane; } sides4_half0 { type symmetryPlane; } sides5_half0 { type symmetryPlane; } out1_half0 { type inletOutlet; inletValue uniform (0 0 0); value uniform (0 0 0); } out2_half0 { type inletOutlet; inletValue uniform (0 0 0); value uniform (0 0 0); } out3_half0 { type inletOutlet; inletValue uniform (0 0 0); value uniform (0 0 0); } out4_half0 { type inletOutlet; inletValue uniform (0 0 0); value uniform (0 0 0); } out5_half0 { type inletOutlet; inletValue uniform (0 0 0); value uniform (0 0 0); } sides5_half1 { type symmetryPlane; } sides4_half1 { type symmetryPlane; } sides3_half1 { type symmetryPlane; } sides2_half1 { type symmetryPlane; } sides1_half1 { type symmetryPlane; } in1_half1 { type turbulentInlet; referenceField uniform (9 0 0); fluctuationScale (0.02 0.01 0.01); value uniform (9 0 0); } in2_half1 { type turbulentInlet; referenceField uniform (9 0 0); fluctuationScale (0.02 0.01 0.01); value uniform (9 0 0); } in3_half1 { type turbulentInlet; referenceField uniform (9 0 0); fluctuationScale (0.02 0.01 0.01); value uniform (9 0 0); } in4_half1 { type turbulentInlet; referenceField uniform (9 0 0); fluctuationScale (0.02 0.01 0.01); value uniform (9 0 0); } in5_half1 { type turbulentInlet; referenceField uniform (9 0 0); fluctuationScale (0.02 0.01 0.01); value uniform (9 0 0); } building { type fixedValue; value uniform ( 0 0 0 ); } } // ************************************************************************* ///*--------------------------------*- C++ -*----------------------------------*\ | ========= | | | \\ / F ield | OpenFOAM: The Open Source CFD Toolbox | | \\ / O peration | Version: 2.4.0 | | \\ / A nd | Web: www.OpenFOAM.org | | \\/ M anipulation | | \*---------------------------------------------------------------------------*/ FoamFile { version 2.0; format ascii; class volScalarField; location "0"; object nuSgs; } // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // dimensions [ 0 2 -1 0 0 0 0 ]; internalField uniform 0; boundaryField { bottomWall { type nutUSpaldingWallFunction; value uniform 0.0; } topWall { type symmetryPlane; } sides1_half0 { type symmetryPlane; } sides2_half0 { type symmetryPlane; } sides3_half0 { type symmetryPlane; } sides4_half0 { type symmetryPlane; } sides5_half0 { type symmetryPlane; } out1_half0 { type inletOutlet; inletValue uniform 0; value uniform 0; } out2_half0 { type inletOutlet; inletValue uniform 0; value uniform 0; } out3_half0 { type inletOutlet; inletValue uniform 0; value uniform 0; } out4_half0 { type inletOutlet; inletValue uniform 0; value uniform 0; } out5_half0 { type inletOutlet; inletValue uniform 0; value uniform 0; } sides5_half1 { type symmetryPlane; } sides4_half1 { type symmetryPlane; } sides3_half1 { type symmetryPlane; } sides2_half1 { type symmetryPlane; } sides1_half1 { type symmetryPlane; } in1_half1 { type fixedValue; value uniform 0.0; } in2_half1 { type fixedValue; value uniform 0.0; } in3_half1 { type fixedValue; value uniform 0.0; } in4_half1 { type fixedValue; value uniform 0.0; } in5_half1 { type fixedValue; value uniform 0.0; } building { type nutUSpaldingWallFunction; value uniform 0.0; } } // ************************************************************************* ///*--------------------------------*- C++ -*----------------------------------*\ | ========= | | | \\ / F ield | OpenFOAM: The Open Source CFD Toolbox | | \\ / O peration | Version: 2.4.0 | | \\ / A nd | Web: www.OpenFOAM.org | | \\/ M anipulation | | \*---------------------------------------------------------------------------*/ FoamFile { version 2.0; format ascii; class volScalarField; location "0"; object k; } // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // dimensions [ 0 2 -2 0 0 0 0 ]; internalField uniform 0; boundaryField { bottomWall { type fixedValue; value uniform 0; } topWall { type symmetryPlane; } sides1_half0 { type symmetryPlane; } sides2_half0 { type symmetryPlane; } sides3_half0 { type symmetryPlane; } sides4_half0 { type symmetryPlane; } sides5_half0 { type symmetryPlane; } out1_half0 { type inletOutlet; inletValue uniform 0; value uniform 0; } out2_half0 { type inletOutlet; inletValue uniform 0; value uniform 0; } out3_half0 { type inletOutlet; inletValue uniform 0; value uniform 0; } out4_half0 { type inletOutlet; inletValue uniform 0; value uniform 0; } out5_half0 { type inletOutlet; inletValue uniform 0; value uniform 0; } sides5_half1 { type symmetryPlane; } sides4_half1 { type symmetryPlane; } sides3_half1 { type symmetryPlane; } sides2_half1 { type symmetryPlane; } sides1_half1 { type symmetryPlane; } in1_half1 { type fixedValue; value uniform 0; } in2_half1 { type fixedValue; value uniform 0; } in3_half1 { type fixedValue; value uniform 0; } in4_half1 { type fixedValue; value uniform 0; } in5_half1 { type fixedValue; value uniform 0; } building { type fixedValue; value uniform 0; } } // ************************************************************************* //```/*--------------------------------*- C++ -*----------------------------------*\ | ========= | | | \\ / F ield | OpenFOAM: The Open Source CFD Toolbox | | \\ / O peration | Version: 2.4.0 | | \\ / A nd | Web: www.OpenFOAM.org | | \\/ M anipulation | | \*---------------------------------------------------------------------------*/ FoamFile { version 2.0; format ascii; class volScalarField; location "0"; object cS; } // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // dimensions [ 0 0 0 0 0 0 0 ]; internalField uniform 0; boundaryField { bottomWall { type zeroGradient; } topWall { type symmetryPlane; } sides1_half0 { type symmetryPlane; } sides2_half0 { type symmetryPlane; } sides3_half0 { type symmetryPlane; } sides4_half0 { type symmetryPlane; } sides5_half0 { type symmetryPlane; } out1_half0 { type zeroGradient; } out2_half0 { type zeroGradient; } out3_half0 { type zeroGradient; } out4_half0 { type zeroGradient; } out5_half0 { type zeroGradient; } sides5_half1 { type symmetryPlane; } sides4_half1 { type symmetryPlane; } sides3_half1 { type symmetryPlane; } sides2_half1 { type symmetryPlane; } sides1_half1 { type symmetryPlane; } in1_half1 { type zeroGradient; } in2_half1 { type zeroGradient; } in3_half1 { type zeroGradient; } in4_half1 { type zeroGradient; } in5_half1 { type zeroGradient; } building { type zeroGradient; } } // ************************************************************************* //solid flange facet normal 0 0 -1 outer loop vertex 1.35 0.00 1.25 vertex 1.35 0.40 1.25 vertex 1.45 0.40 1.25 endloop endfacet facet normal 0 0 -1 outer loop vertex 1.35 0.00 1.25 vertex 1.45 0.40 1.25 vertex 1.45 0.00 1.25 endloop endfacet facet normal -1 0 0 outer loop vertex 1.35 0.00 1.25 vertex 1.35 0.00 1.15 vertex 1.35 0.40 1.15 endloop endfacet facet normal -1 0 0 outer loop vertex 1.35 0.00 1.25 vertex 1.35 0.40 1.15 vertex 1.35 0.40 1.25 endloop endfacet facet normal 0 -1 0 outer loop vertex 1.35 0.00 1.25 vertex 1.45 0.00 1.15 vertex 1.45 0.00 1.25 endloop endfacet facet normal 0 -1 0 outer loop vertex 1.35 0.00 1.25 vertex 1.45 0.00 1.15 vertex 1.35 0.00 1.15 endloop endfacet facet normal 0 0 1 outer loop vertex 1.35 0.00 1.15 vertex 1.45 0.00 1.15 vertex 1.45 0.40 1.15 endloop endfacet facet normal 0 0 1 outer loop vertex 1.35 0.00 1.15 vertex 1.45 0.40 1.15 vertex 1.35 0.40 1.15 endloop endfacet facet normal 1 0 0 outer loop vertex 1.45 0.00 1.25 vertex 1.45 0.40 1.25 vertex 1.45 0.40 1.15 endloop endfacet facet normal 1 0 0 outer loop vertex 1.45 0.00 1.25 vertex 1.45 0.40 1.15 vertex 1.45 0.00 1.15 endloop endfacet facet normal 0 1 0 outer loop vertex 1.35 0.40 1.25 vertex 1.35 0.40 1.15 vertex 1.45 0.40 1.15 endloop endfacet facet normal 0 1 0 outer loop vertex 1.35 0.40 1.25 vertex 1.45 0.40 1.25 vertex 1.45 0.40 1.15 endloop endfacet endsolid flange -
非常感谢 @gtian 提供的模型，最近和一个师弟在模拟TPU建筑时，迎风面风压总是计算整体偏小。我按照提供的模型及边界条件设置再计算了下，有几个疑问，再请教一下：
1. 速度入口目标值选择问题
描述：提供的边界条件选用了
turbulentInlet,参数设定为referenceField uniform (9 0 0); fluctuationScale (0.02 0.01 0.01);，这是一个平均速度为 9m/s 的均匀湍流场。但是TPU实验用的是指数率平均风速剖面，如下图所示。根据TPU数据库，对于建筑模型1:1:4，风速剖面指数为1/4的实验，参考风速为11m/s，即：$U(z)=U_{ref}\left ( \frac{z}{H_{ref}} \right )^{\alpha},U_{ref}=11m/s,H_{ref}=0.4m,\alpha =1/4$
问题：入口速度目标值不一致，计算结果跟TPU的风洞实验对比，是否会不太合适？
2. 平均风压系数的计算问题
描述：平均风压系数计算公式为：$ \bar{C}_{p}= \frac{\bar{p}}{\frac{1}{2} \rho U_H^2} $，其中 $\bar{p}$ 为平均压力，$\rho$为空气密度，$U_H$ 为建筑物高度处风速。对于TPU指数率剖面，本工况取11m/s。
问题：对于 9m/s 的均匀湍流场，$U_H$是如何确定？以及由OpenFOAM导出的压力值，计算风压系数时的分母是否需要除以密度$\rho$？3. 计算网格边界层设置问题
在计算钝体绕流问题时，计算域的底面、建筑物表面的边界层，是否需设置到符合壁面函数的
y+大小？如果没设置边界层，用了nutUSpaldingWallFunction壁面函数，对结果影响有多大？4. 根据 @gtian 提供的边界条件及模型，我们的计算结果
我们采用
Smagorinsky模型计算，其他的边界条件及求解格式都与提供算例一致，也用的是 9m/s 的均匀湍流场。从1s后开始统计平均场，总共计算24s。计算结果如下。可以看出最大平均风压为41Pa，计算最大风压系数约为：$41/(0.5\times 1.225\times 9\times 9)=0.826$。
5. 前期我们模拟TPU114模型的工况设置及结果
前期我们模拟TPU114模型设置，见下载地址：
链接：https://pan.baidu.com/s/161fNMd9disZoZtW-CZgaRg
提取码：rznq用fluentmeshing画的网格，采用groovyBC加载平均风剖面计算，$U(z)=U_{ref}\left ( \frac{z}{H_{ref}} \right )^{\alpha},U_{ref}=11m/s,H_{ref}=0.4m,\alpha =1/4$。
（1）建筑物平均风压结果：最大平均风压系数为 $54/(1/2\times 1.225\times 11\times 11)=0.728$，而实验值的最大平均风压系数为0.874，差距0.15。
（2）平均速度剖面：
（3）OpenFOAM与Fluent结果对比
①迎风面平均风压云图
②2/3H高度处的平均风压结果。OpenFOAM迎风面平均压力系数结果偏小0.1几，但fluent结果与实验值吻合程度非常好。其中两个侧面和背风面，有旋涡脱落，没有加脉动风，因此模拟结果较差是正常的。由于迎风面结果与平均风剖面大小是比较相关的，正常只加平均风，应该模拟结果是最好的。但目前OpenFOAM结果较差。
③除提供的算例设置外，我们尝试过修改网格、边界条件、湍流模型、fvSchemes、fvSolution，但都是这样的结果。可能一直认为加载平均风，OpenFOAM的迎风面平均风压系数应该和Fluent一样很吻合才比较合理，但一直没找到原因。
本来应该是一个很成熟的计算算例，但实在不知道什么因素影响了结果，请各位大佬指点迷津，实在感谢
-
@coolhhh 我之前没有算过风压这类问题，还挺有兴趣的。我对这个u_h，也是非常疑惑。我用的u_h是柱体中心横截面上（u_h）_max。之所以有疑惑，是关于这个u_h的选择位置的，因为我们一般算cp用到是远场处的速度（我之前做过3900圆柱绕流，https://doi.org/10.1063/5.0012358）。事实上，柱体上方会有一个强剪切区，流动分离后会形成回流区。柱体侧边同理。因此，我很疑惑的是u_h 是如何定义的？ 是 h高度处的最大速度吗？如果有流动分离，即使在初始边界上设置了u_h，那么这个u_h也会发生变化啊。反倒是都用，2/3高度处的流场数据（pmean-0）/ （0.5* u_0.667 ^2），我觉得更合理。因为不可压流p就是p/rho，可以看一下文件的量纲，所以就不需要在除密度了。
-
@coolhhh 还有，尽量不要用smagorinsky做亚格子模型，因为这个亚格子模型中的系数C是个定值。他会有两个问题，1.边界上，会overestimate边界附近的粘度。2. 是全场是一个固定系数的话，会忽视掉能量的反级串backscatter现象。one-equation或者dynamic的亚格子模型都比smagorinsky要好一些。对于静压力，我试一下你给的那个边界条件跑一下看看。还有就是关于总用时的。我一般分成两部分，一部分是从初识态到流场稳定（200-300个t*）。第二阶段，统计过程，100-300个t*。t*=特征长度/特征速度。比如H/u*，摩擦速度可以用sqrt（压差阻力/横面积）得到。
