钝体建筑扰流的大涡模拟
-
\*---------------------------------------------------------------------------*/ FoamFile { version 2.0; format ascii; class volVectorField; location "1"; object U; } // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // dimensions [0 1 -1 0 0 0 0]; internalField uniform (0 0 0); boundaryField { inlet { type timeVaryingMappedFixedValue; offset (0 0 0); setAverage off; } outlet { type zeroGradient; } upperwall { type symmetry; } bottom { type noSlip; } symmetry { type symmetry; } building { type noSlip; } }// ************************************************************************* // /*--------------------------------*- C++ -*----------------------------------*\ | ========= | | | \\ / F ield | OpenFOAM: The Open Source CFD Toolbox | | \\ / O peration | Version: v2012 | | \\ / A nd | Website: www.openfoam.com | | \\/ M anipulation | | \*---------------------------------------------------------------------------*/ FoamFile { version 2.0; format ascii; class volScalarField; location "1"; object p; } // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // dimensions [0 2 -2 0 0 0 0]; internalField uniform 0; boundaryField { inlet { type zeroGradient; } outlet { type fixedValue; value uniform 0; //type zeroGradient; } upperwall { type symmetry; } bottom { type zeroGradient; } symmetry { type symmetry; } building { type zeroGradient; } }// ************************************************************************* // /*--------------------------------*- C++ -*----------------------------------*\ | ========= | | | \\ / F ield | OpenFOAM: The Open Source CFD Toolbox | | \\ / O peration | Version: v2012 | | \\ / A nd | Website: www.openfoam.com | | \\/ M anipulation | | \*---------------------------------------------------------------------------*/ FoamFile { version 2.0; format ascii; class volScalarField; location "1"; object nut; } // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // dimensions [0 2 -1 0 0 0 0]; internalField uniform 0; boundaryField { inlet { type calculated; value uniform 0; } outlet { type calculated; value uniform 0; } upperwall { type symmetry; } bottom { type nutUSpaldingWallFunction; value $internalField; //type fixedValue; //value uniform 0; } symmetry { type symmetry; } building { type nutUSpaldingWallFunction; value $internalField; //type fixedValue; //value uniform 0; } }// ************************************************************************* // /*--------------------------------*- C++ -*----------------------------------*\ | ========= | | | \\ / F ield | OpenFOAM: The Open Source CFD Toolbox | | \\ / O peration | Version: v2012 | | \\ / A nd | Website: www.openfoam.com | | \\/ M anipulation | | \*---------------------------------------------------------------------------*/ FoamFile { version 2.0; format ascii; class dictionary; location "constant"; object turbulenceProperties; } // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // simulationType LES; LES { LESModel WALE; turbulence on; printCoeffs on; delta cubeRootVol; cubeRootVolCoeffs { deltaCoeff 1; } PrandtlCoeffs { delta cubeRootVol; cubeRootVolCoeffs { deltaCoeff 1; } smoothCoeffs { delta cubeRootVol; cubeRootVolCoeffs { deltaCoeff 1; } maxDeltaRatio 1.1; } Cdelta 0.158; } vanDriestCoeffs { delta cubeRootVol; cubeRootVolCoeffs { deltaCoeff 1; } smoothCoeffs { delta cubeRootVol; cubeRootVolCoeffs { deltaCoeff 1; } maxDeltaRatio 1.1; } Aplus 26; Cdelta 0.158; } smoothCoeffs { delta cubeRootVol; cubeRootVolCoeffs { deltaCoeff 1; } maxDeltaRatio 1.1; } }// ************************************************************************* // /*--------------------------------*- C++ -*----------------------------------*\ | ========= | | | \\ / F ield | OpenFOAM: The Open Source CFD Toolbox | | \\ / O peration | Version: 4.x | | \\ / A nd | Web: www.OpenFOAM.org | | \\/ M anipulation | | \*---------------------------------------------------------------------------*/ FoamFile { version 2.0; format ascii; class dictionary; location "system"; object fvSchemes; } // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // //shijianyijiedushuxiang ddtSchemes // { default backward; } //tidu gradSchemes { default Gauss linear; //grad(nuTilda) cellLimited Gauss linear 1; //grad(U) cellLimited Gauss linear 1; //grad(P) Gauss linear; } //sandduxiang divSchemes ///////////////////////////////////////////////////// { default none; div(phi,U) bounded Gauss linearUpwind grad(U); div(phi,k) bounded Gauss linearUpwind grad(k); div(phi,B) bounded Gauss linearUpwind grad(B); div(phi,omega) bounded Gauss linearUpwind grad(omega); div(phi,B) bounded Gauss linear; //div(phi,nuTilda) Gauss limitedLinear 1; div((nuEff*dev2(T(grad(U))))) Gauss linear; } laplacianSchemes // { default Gauss linear limited corrected 0.33; } interpolationSchemes // { default linear; } snGradSchemes // { default limited corrected 0.33; } wallDist { method meshWave; }// ************************************************************************* // /*--------------------------------*- C++ -*----------------------------------*\ | ========= | | | \\ / F ield | OpenFOAM: The Open Source CFD Toolbox | | \\ / O peration | Version: 4.x | | \\ / A nd | Web: www.OpenFOAM.org | | \\/ M anipulation | | \*---------------------------------------------------------------------------*/ FoamFile { version 2.0; format ascii; class dictionary; location "system"; object fvSolution; } // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // solvers { p { solver GAMG; tolerance 1e-6; relTol 0.1; smoother GaussSeidel; //nCellsInCoarsestLevel 50; } pFinal { $p; tolerance 1e-6; relTol 0; } U { solver PBiCG; preconditioner DILU; tolerance 1e-6; relTol 0.05; } UFinal { solver PBiCG; preconditioner DILU; tolerance 1e-6; relTol 0; } } PISO { nCorrectors 3; nNonOrthogonalCorrectors 1; //pRefCell 0; pRefPoint (2 0 1.5); pRefValue 0; } relaxationFactors { fields { p 0.3; } equations { U 1; } //"U.*" 0.5; //"nuTilda.*" 1; } // ************************************************************************* //李老师这是我的设置，感谢各位老师朋友的解答
-
snappyHexMeshDict 中把snap 和 add 关掉，就会生成结构化的六面体网格，网格质量会号很多， 只需要有stl文件即可.
// Which of the steps to run castellatedMesh true; snap false; addLayers false; geometry { flange11.stl { type triSurfaceMesh; name flange11; }stl文件如下：
solid flange facet normal 0 0 -1 outer loop vertex 0.170 0.000 0.230 vertex 0.170 0.020 0.230 vertex 0.190 0.020 0.230 endloop endfacet facet normal 0 0 -1 outer loop vertex 0.170 0.000 0.230 vertex 0.190 0.020 0.230 vertex 0.190 0.000 0.230 endloop endfacet facet normal -1 0 0 outer loop vertex 0.170 0.000 0.230 vertex 0.170 0.000 0.210 vertex 0.170 0.020 0.210 endloop endfacet facet normal -1 0 0 outer loop vertex 0.170 0.000 0.230 vertex 0.170 0.020 0.210 vertex 0.170 0.020 0.230 endloop endfacet facet normal 0 -1 0 outer loop vertex 0.170 0.000 0.230 vertex 0.190 0.000 0.210 vertex 0.190 0.000 0.230 endloop endfacet facet normal 0 -1 0 outer loop vertex 0.170 0.000 0.230 vertex 0.190 0.000 0.210 vertex 0.170 0.000 0.210 endloop endfacet facet normal 0 0 1 outer loop vertex 0.170 0.000 0.210 vertex 0.190 0.000 0.210 vertex 0.190 0.020 0.210 endloop endfacet facet normal 0 0 1 outer loop vertex 0.170 0.000 0.210 vertex 0.190 0.020 0.210 vertex 0.170 0.020 0.210 endloop endfacet facet normal 1 0 0 outer loop vertex 0.190 0.000 0.230 vertex 0.190 0.020 0.230 vertex 0.190 0.020 0.210 endloop endfacet facet normal 1 0 0 outer loop vertex 0.190 0.000 0.230 vertex 0.190 0.020 0.210 vertex 0.190 0.000 0.210 endloop endfacet facet normal 0 1 0 outer loop vertex 0.170 0.020 0.230 vertex 0.170 0.020 0.210 vertex 0.190 0.020 0.210 endloop endfacet facet normal 0 1 0 outer loop vertex 0.170 0.020 0.230 vertex 0.190 0.020 0.230 vertex 0.190 0.020 0.210 endloop endfacet endsolid flange -
@hitsc30 如果你提到的边界层指的是wall周围的对数区以下的粘性子层和buffer layer的话。可以有两种方式，一种是wall model LES，通过壁面模型来调nut，这时第一层网格在y+=30 即可。另一种是wall resolved LES，这种情况只能加密网格。blockmesh和shm都可以加密网格，blockmesh可以控制aspect ratio来加密，shm是4等分加密。这种情况下如果要resolve到粘性子层（viscous sublayer），需要在y+<5 以内至少有三个点，也就是第一点y+<1.
-
@hitsc30 你这个涉及的东西有点多，排查起来会比较复杂。首先建议抛开fluent里的设置，因为大涡模拟这块openfoam和fluent差的还有点多。
首先是网格，snappy如果不好生成的话，试试cfmesh，这个比snappy好用，比fluentmeshing我感觉也就是边界层生成差一些。考虑到你这个几何相对应该比较规整，cfmesh应该可以的。
其次是入口条件。openfoam大涡模拟的入口条件这块是不如fluent的vortex method的，fluent的很好用，不怎么发散，生成的湍流也很强。我看你用的是timeVaryingMappedFixedValue，这是用的实验数据或者测试数据吗？你fluent是怎么用的？
然后是nut的边界条件，我之前也是发现nutUSpaldingWallFunction很好用，但我不确定它应该和WALE一起。WALE一般是传热领域用的多，原因是这个模型是不需要壁面函数的，而传热主要就是要解析壁面，不能模化。WALE感觉和你加的nut边界条件可能会冲突？
再然后是数值格式，一个是pisoFoam，openfoam里面的pisoFoam外迭代只有一步，可能收敛不太好。试试pimpleFoam加上外迭代吧，这个其实你在fluent里面应该也设置过。或者把你目前计算的残差曲线贴出来看看。
最后，为了排查这些问题，建议先跑个非稳态的RANS试试，不然理不出头绪来。其实最好先跑个稳态RANS，没问题了再非稳态RANS，再没问题了再上LES。
-
@gtian 在 钝体建筑扰流的大涡模拟 中说：
@hitsc30 如果你提到的边界层指的是wall周围的对数区以下的粘性子层和buffer layer的话。可以有两种方式，一种是wall model LES，通过壁面模型来调nut，这时第一层网格在y+=30 即可。另一种是wall resolved LES，这种情况只能加密网格。blockmesh和shm都可以加密网格，blockmesh可以控制aspect ratio来加密，shm是4等分加密。这种情况下如果要resolve到粘性子层（viscous sublayer），需要在y+<5 以内至少有三个点，也就是第一点y+<1.
感谢回复，这个确实可以在后续重点突击一下。此前我曾用blockmesh做网格背景、snappyhexmesh加密做过一个网格。当时我得到的计算值也和我用fluentmeshing商软所做的商业网格差不太多，没太考虑网格的质量。您提出的这个是我当前需要考虑的一个方向。我当前觉得目前比较影响我接过的可能是离散格式，在这方面的知识储备非常的欠缺，只能选一些高阶离散格式，可能对结果有不小的影响。 -
@cccrrryyy 您好，感谢回复。目前接触到的在openfoam中做建筑风荷载模拟的论文主要是tamura教授在北京交通大学和重庆大学领导的学科组推出来的，都采用的是pisofoam算法。湍流风入口是由CDRFG方法生成的，在中文网有个帖子谈论过这种植入。也咨询过北交的师兄，他们植入的CDRFG算法也是类似的植入，该方法应该是没有问题的。
nut函数，我也和您的观点一样。该函数能自动切换成粘性层的线性表达式和对数表达式。再者是您提出的WALE，这个确实是需要改动一下.目前接触的论文中，foam中都是用的经典的smagorinksy。fluent中，国内硕士很多用的是WALE结果对的确实不是很好。您提到的WALE这个模型的具体作用，对代码认识比较浅，确实之前没考虑过，确实是我需要再次调整的方向。
当前的两个重点是离散格式和亚格子模型，然后是上面的师兄说的网格的质量。
ps：感觉foam计算y+很奇怪，我的边界层网格每次改变的不少，感觉y+变化趋势和我想的南辕北辙。
-
@hitsc30
y+ 的理解不复杂，看code需要一点点的耐心就可以。我给你贴一下yPlusLES。if (isA<wallFvPatch>(currPatch)) { yPlus.boundaryField()[patchi] = d[patchi] *sqrt ( nuEff.boundaryField()[patchi] *mag(U.boundaryField()[patchi].snGrad()) ) /nuLam.boundaryField()[patchi]; const scalarField& Yp = yPlus.boundaryField()[patchi]; Info<< "Patch " << patchi << " named " << currPatch.name() << " y+ : min: " << gMin(Yp) << " max: " << gMax(Yp) << " average: " << gAverage(Yp) << nl << endl; }
1)
2)
3)我写的，你可以输出来看看每一个量的大小
if (isA<wallFvPatch>(currPatch)) { frictionvelocityLES.boundaryField()[patchi] = sqrt ( nuEff.boundaryField()[patchi] *mag(U.boundaryField()[patchi].snGrad()) ); const scalarField& Yp = frictionvelocityLES.boundaryField()[patchi]; Info<< "Patch " << patchi << " named " << currPatch.name() << " friction velocity u* : min: " << gMin(Yp) << " max: " << gMax(Yp) << " average: " << gAverage(Yp) << nl << endl; Info<< "Patch " << patchi << " named " << currPatch.name() << "nuLam : min: " << gMin(nuLam) << " max: " << gMax(nuLam) << " average: " << gAverage(nuLam) << nl << endl; Info<< "Patch " << patchi << " named " << currPatch.name() << "d : min: " << gMin(d) << " max: " << gMax(d) << " average: " << gAverage(d) << nl << endl; }近壁面第一层网格的话 ，d= 0.5* 第一个cell的长度. du/dy = （u_infty - 0）/ (d-0)
nuefff = nu （层流 laminar） + nut（wall 上 nut =0 或者是你给的壁面模型的值）我给的nut在壁面上就是0.
object nuSgs; } // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // dimensions [ 0 2 -1 0 0 0 0 ]; internalField uniform 0; boundaryField { bottomWall { type fixedValue; value uniform 0; } -
@gtian 在 钝体建筑扰流的大涡模拟 中说：
@hitsc30
y+ 的理解不复杂，看code需要一点点的耐心就可以。我给你贴一下yPlusLES。if (isA<wallFvPatch>(currPatch)) { yPlus.boundaryField()[patchi] = d[patchi] *sqrt ( nuEff.boundaryField()[patchi] *mag(U.boundaryField()[patchi].snGrad()) ) /nuLam.boundaryField()[patchi]; const scalarField& Yp = yPlus.boundaryField()[patchi]; Info<< "Patch " << patchi << " named " << currPatch.name() << " y+ : min: " << gMin(Yp) << " max: " << gMax(Yp) << " average: " << gAverage(Yp) << nl << endl; }
我写的，你可以输出来看看每一个量的大小
if (isA<wallFvPatch>(currPatch)) { frictionvelocityLES.boundaryField()[patchi] = sqrt ( nuEff.boundaryField()[patchi] *mag(U.boundaryField()[patchi].snGrad()) ); const scalarField& Yp = frictionvelocityLES.boundaryField()[patchi]; Info<< "Patch " << patchi << " named " << currPatch.name() << " friction velocity u* : min: " << gMin(Yp) << " max: " << gMax(Yp) << " average: " << gAverage(Yp) << nl << endl; Info<< "Patch " << patchi << " named " << currPatch.name() << "nuLam : min: " << gMin(nuLam) << " max: " << gMax(nuLam) << " average: " << gAverage(nuLam) << nl << endl; Info<< "Patch " << patchi << " named " << currPatch.name() << "d : min: " << gMin(d) << " max: " << gMax(d) << " average: " << gAverage(d) << nl << endl; }近壁面第一层网格的话 ，d= 0.5* 第一个cell的长度. du/dy = （u_infty - 0）/ (d-0)
nuefff = nu （层流 laminar） + nut（wall 上 nut =0 或者是你给的壁面模型的值）我给的nut在壁面上就是0.
object nuSgs; } // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // dimensions [ 0 2 -1 0 0 0 0 ]; internalField uniform 0; boundaryField { bottomWall { type fixedValue; value uniform 0; }非常感谢！
-
-
@李东岳 李老师我上点图给您和各位老师看一下。我们使用了fluent和openfoam进行平均风剖面的模拟。下图是我们得到的结果
下述是平均风剖面对比:
而后是我用脉动算出来的结果，分别是平均风压和脉动风压
下述是我的实验数据迎风面风压和模拟的迎风面风压
误差有点离谱
当前是按照tamura教授的设置帮边界条件、离散格式、求解器、亚格子模型全换了。就是采用的网格和他们的有所区别（他们采用的snappyhexmesh），准备搞个和他们一样的网格试一试。实在不行就只能放弃了
主要现在用foam用了半年，从一个啥都不知道的小白在师兄和论坛的各位前辈带着下，慢慢的进步对这个东西用习惯了，人菜瘾大。就想证实到底是我设置什么的不行，还是foam可能这方面存在改进。就是最后放弃foam改用fluent算也不让自己后悔。
抱歉，李博。这几天电脑卡的要死，开启fluentmeshing就卡死，一直想着要回复您。但是也没法打开软件。这是我再fluentmeshing中检测的问题。我当前的问题主要是，我画过粗略的多面体和精细的多面体网格是否对于结果优化都没太大的影响了。当前想着试一试snappyhexmesh画六面体拯救一下我，在使用中发现其对于边界层网格生成有点随缘化，网格在90度拐角处处理的不是很好。
