本案例利用OpenFOAM中的simpleFoam求解器计算建筑物外流场。
注:案例位于incompressiblesimpleFoamwind
”
1 计算模型
建筑物模型如下图所示。
外流计算域在blockMeshDict
中进行指定。外流场几何尺寸x轴方向[-20, 330],y方向[-50,230],z方向[0,140]。
本算例计算网格采用snappyHexMesh进行生成。其实也可以利用第三方网格工具生成后导入到OpenFOAM中。
2 计算网格
网格生成命令为:
surfaceFeatures
blockMesh
snappyHexMesh -overwrite
网格生成完毕后如下图所示。
2.1 构造面特征
在进行网格生成之前,需要先利用程序surfaceFeatures
构造几何特征。该程序需要利用字典文件surfaceFeaturesDict
,此字典文件位于system
路径下,文件内容为:
FoamFile
{
version 2.0;
format ascii;
class dictionary;
object surfaceFeaturesDict;
}
// 这里指定算例几何文件
surfaces ("buildings.obj");
// 添加面特征提取信息
#includeEtc "caseDicts/surface/surfaceFeaturesDict.cfg"
注:在实际应用过程中,可以使用命令
foamGet surfaceFeatures
直接添加字典文件,之后修改字典文件中的算例几何信息即可。”
字典文件中包含的surfaceFeaturesDict.cfg
文件位于$WM_PROJECT_DIR/etc/caseDicts/surface
文件夹中。文件内容如下所示,此文件一般情况下不需要修改。
includedAngle 150;
subsetFeatures
{
nonManifoldEdges yes;
openEdges yes;
}
trimFeatures
{
minElem 0;
minLen 0;
}
writeObj yes;
运行surfaceFeatures
命令后会在constant/triSurface
文件夹下创建一个名为buildings.eMesh
的文件,该文件中包含了所有释放出的几何面信息,在后续的snappyHexMesh
生成网格时需要使用到此文件。
2.2 生成背景网格
背景网格使用blockMesh
生成。
blockMesh命令需要配合blockMeshDict
字典文件一起运行。本算例的blockMeshDict
字典文件内容如下所示:
FoamFile
{
version 2.0;
format ascii;
class dictionary;
object blockMeshDict;
}
// * * * * * * * * * * * * * * * * * * * * //
// 定义计算域尺寸及网格生成所需的数据
backgroundMesh
{
xMin -20; // L = 350
xMax 330;
yMin -50; // L = 280
yMax 230;
zMin 0;
zMax 140;
xCells 25;
yCells 20;
zCells 10;
}
convertToMeters 1;
// 本算例外流场是一个长方体
// 这里指定长方体的六个定点坐标
vertices
(
($:backgroundMesh.xMin $:backgroundMesh.yMin $:backgroundMesh.zMin)
($:backgroundMesh.xMax $:backgroundMesh.yMin $:backgroundMesh.zMin)
($:backgroundMesh.xMax $:backgroundMesh.yMax $:backgroundMesh.zMin)
($:backgroundMesh.xMin $:backgroundMesh.yMax $:backgroundMesh.zMin)
($:backgroundMesh.xMin $:backgroundMesh.yMin $:backgroundMesh.zMax)
($:backgroundMesh.xMax $:backgroundMesh.yMin $:backgroundMesh.zMax)
($:backgroundMesh.xMax $:backgroundMesh.yMax $:backgroundMesh.zMax)
($:backgroundMesh.xMin $:backgroundMesh.yMax $:backgroundMesh.zMax)
);
// 创建一个block,并指定沿三个方向的网格数量
blocks
(
hex (0 1 2 3 4 5 6 7)
(
$:backgroundMesh.xCells
$:backgroundMesh.yCells
$:backgroundMesh.zCells
)
simpleGrading (1 1 1)
);
edges
(
);
// 下面定义边界
// 创建了4个边界:inletoutletgroundfrontAndBack
boundary
(
inlet
{
type patch;
faces
(
(0 3 7 4)
);
}
outlet
{
type patch;
faces
(
(1 5 6 2)
);
}
ground
{
type wall;
faces
(
(0 1 2 3)
);
}
// 边界类型为symmetry
frontAndBack
{
type symmetry;
faces
(
(0 4 5 1)
(3 2 6 7)
(4 7 6 5)
);
}
);
mergePatchPairs
(
);
2.3 生成体网格
体网格利用snappyHexMesh
生成,需要在system
文件夹中定义字典文件snappyHexMeshDict
。本算例中该文件定义如下所示。
FoamFile
{
version 2.0;
format ascii;
class dictionary;
object snappyHexMeshDict;
}
// * * * * * * * * * * * * * * * * * * * //
// 配置文件,一般不改
#includeEtc "caseDicts/mesh/generation/snappyHexMeshDict.cfg"
castellatedMesh on;
snap on;
addLayers off;
geometry
{
buildings
{
type triSurfaceMesh;
file "buildings.obj"; //指定几何文件
}
refinementBox
{
type searchableBox;
min ( 0 0 0);
max (250 180 90);
}
};
castellatedMeshControls
{
features
(
{ file "buildings.eMesh"; level 1; } //指定几何面信息文件
);
refinementSurfaces
{
buildings
{
level (3 3);
patchInfo { type wall; }
}
}
refinementRegions
{
refinementBox
{
mode inside;
levels ((1E15 2));
}
}
locationInMesh (1 1 1);
}
snapControls
{
explicitFeatureSnap true;
implicitFeatureSnap false;
}
addLayersControls
{
layers
{
"CAD.*"
{
nSurfaceLayers 2;
}
}
relativeSizes true;
expansionRatio 1.2;
finalLayerThickness 0.5;
minThickness 1e-3;
}
meshQualityControls
{}
writeFlags
(
// scalarLevels
// layerSets
// layerFields
);
mergeTolerance 1e-6;
snappyHexMeshDict文件较为复杂,该文件模板位于$WM_PROJECT_DIR/etc/caseDicts/mesh/generation
文件中,应用过程中可以使用命令foamGet snappyHexMeshDict
将该文件直接添加到system文件夹中,之后修改文件中与算例几何相关的文件名即可。
3 物理模型
本算例计算的是不可压缩湍流流动,未考虑传热。因此仅需在constant
文件夹中定义属性文件transportProperties
及momentumTransport
即可。
注:在早期OpenFOAM版本中,momentumTransport文件名为turbulenceProperties。
”
-
在 transportProperties
文件中定义流体物性,如下所示:
FoamFile
{
version 2.0;
format ascii;
class dictionary;
object transportProperties;
}
// * * * * * * * * * * * * * * * * * * * //
transportModel Newtonian;
// 定义运动粘度为1.5e-5 m2/s
nu [0 2 -1 0 0 0 0] 1.5e-05;
-
在 momentumTransport
文件中定义湍流模型。本算例采用kEpsilon
湍流模型。文件内容如下:
FoamFile
{
version 2.0;
format ascii;
class dictionary;
object RASProperties;
}
// * * * * * * * * * * * * * * * * * * * * * //
simulationType RAS;
RAS
{
RASModel kEpsilon;
turbulence on;
printCoeffs on;
}
4 边界条件
边界条件的指定位于0
文件夹,其中需要指定p、U、k、espilon、nut
等物理量。边界名称在文件constant/polyMesh/boundary
中指定,该文件如下所示。
FoamFile
{
version 2.0;
format ascii;
class polyBoundaryMesh;
location "constant/polyMesh";
object boundary;
}
// * * * * * * * * * * * * * * * //
5
(
inlet
{
type patch;
nFaces 494;
startFace 558513;
}
outlet
{
type patch;
nFaces 200;
startFace 559007;
}
ground
{
type wall;
inGroups List 1(wall);
nFaces 6310;
startFace 559207;
}
frontAndBack
{
type symmetry;
inGroups List 1(symmetry);
nFaces 1000;
startFace 565517;
}
buildings
{
type wall;
inGroups List 1(wall);
nFaces 22723;
startFace 566517;
}
)
可以看到,一共定义了5个边界:inlet、outlet、ground、frontAndBack以及buildings。其中边界inlet与outlet类型为patch;frontAndBack的类型为symmetry;ground与buildings边界类型为wall。在边界文件中需要指定这些边界条件。OpenFOAM提供了一些便利条件,可以给相同类型的边界指定统一的边界条件,如wall、symmetry类型的边界。
4.1 p文件
p文件中指定各边界上的压力分布。
FoamFile
{
version 2.0;
format ascii;
class volScalarField;
object p;
}
// * * * * * * * * * * * * * * * * //
dimensions [0 2 -2 0 0 0 0];
internalField uniform 0;
boundaryField
{
inlet
{
type zeroGradient;
}
// 指定边界outlet类型为总压出口
// 不可压缩流体中仅需要指定p0即可
outlet
{
type totalPressure;
p0 uniform 0;
}
// 这里是简化用法
// 等效于"(ground|building)"
wall
{
type zeroGradient;
}
#includeEtc "caseDicts/setConstraintTypes"
}
这个边界条件指定文件中并没有特别复杂的地方,下面两点需要注意。
1、以边界类型指定边界条件
如文件中的:
wall
{
type zeroGradient;
}
其意义为指定所有wall
类型的边界(本算例为ground与buildings)条件。其等效于:
"(ground|buildings)"
{
type zeroGradient;
}
2、包含文件setConstraintTypes
在边界文件的最后一行#includeEtc "caseDicts/setConstraintTypes"
是包含了一个文件setContraintType,该文件位于文件路径/opt/openfoam8/etc/caseDicts/setConstraintTypes
,在该文件中为一些常用的约束边界指定了类型。在本算例中,该语句等效于:
symmetry
{
type symmetry;
}
或写成更直观的形式:
frontAndBack
{
type symmetry;
}
因此p
文件可以改写为:
FoamFile
{
version 2.0;
format ascii;
class volScalarField;
object p;
}
// * * * * * * * * * * * * //
dimensions [0 2 -2 0 0 0 0];
internalField uniform 0;
boundaryField
{
inlet
{
type zeroGradient;
}
outlet
{
type totalPressure;
p0 uniform 0;
}
"(ground|buildings)"
{
type zeroGradient;
}
frontAndBack
{
type symmetry;
}
}
4.3 U文件
U
文件内容如下所示。
FoamFile
{
version 2.0;
format ascii;
class volVectorField;
object U;
}
// * * * * * * * * * * * * * * //
Uinlet (10 0 0);
dimensions [0 1 -1 0 0 0 0];
internalField uniform (0 0 0);
boundaryField
{
// 指定入口速度为X方向10 m/s
inlet
{
type fixedValue;
value uniform $Uinlet;
}
// 指定出口速度类型为pressureInletOutletVelocity
/* 此边界条件适用于指定了压力的压力边界。
对于流出,应用零梯度条件;对于流入,速度由边界面的法向分量获得。
*/
outlet
{
type pressureInletOutletVelocity;
value uniform (0 0 0);
}
wall
{
type noSlip;
}
#includeEtc "caseDicts/setConstraintTypes"
}
4.3 k文件
k文件内容如下所示。
FoamFile
{
version 2.0;
format ascii;
class volScalarField;
object k;
}
// approx k = 1.5*(I*U)^2 ; I = 0.1
kInlet 1.5;
dimensions [0 2 -2 0 0 0 0];
internalField uniform $kInlet;
boundaryField
{
inlet
{
type fixedValue;
value uniform $kInlet;
}
outlet
{
type inletOutlet;
inletValue uniform $kInlet;
value uniform $kInlet;
}
wall
{
type kqRWallFunction;
value uniform $kInlet;
}
#includeEtc "caseDicts/setConstraintTypes"
}
4.4 epsilon文件
epsilon文件如下所示。
FoamFile
{
version 2.0;
format ascii;
class volScalarField;
object epsilon;
}
// Cmu^0.75 * k^1.5 / L ; L =10
epsilonInlet 0.03;
dimensions [0 2 -3 0 0 0 0];
internalField uniform $epsilonInlet;
boundaryField
{
inlet
{
type fixedValue;
value uniform $epsilonInlet;
}
outlet
{
type inletOutlet;
inletValue uniform $epsilonInlet;
value uniform $epsilonInlet;
}
wall
{
type epsilonWallFunction;
value uniform $epsilonInlet;
}
#includeEtc "caseDicts/setConstraintTypes"
}
4.5 nut文件
nut文件内容如下所示。
FoamFile
{
version 2.0;
format ascii;
class volScalarField;
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;
}
wall
{
type nutkWallFunction;
value uniform 0;
}
#includeEtc "caseDicts/setConstraintTypes"
}
5 求解控制
5.1 controlDict文件
controlDict文件如下所示。
FoamFile
{
version 2.0;
format ascii;
class dictionary;
object controlDict;
}
// * * * * * * * * * * * * * //
application simpleFoam;
startFrom latestTime;
startTime 0;
stopAt endTime;
// 迭代计算400步
endTime 400;
deltaT 1;
writeControl timeStep;
writeInterval 50;
purgeWrite 0;
writeFormat ascii;
writePrecision 6;
writeCompression off;
timeFormat general;
timePrecision 6;
runTimeModifiable true;
functions
{
}
5.2 fvSchemes文件
fvSchemes文件指定求解算法。本算例的fvSchemes
文件如下所示。
FoamFile
{
version 2.0;
format ascii;
class dictionary;
object fvSchemes;
}
// 指定采用稳态计算
ddtSchemes
{
default steadyState;
}
gradSchemes
{
default Gauss linear;
limited cellLimited Gauss linear 1;
grad(U) $limited;
grad(k) $limited;
grad(epsilon) $limited;
}
divSchemes
{
default none;
div(phi,U) bounded Gauss linearUpwind limited;
turbulence bounded Gauss limitedLinear 1;
div(phi,k) $turbulence;
div(phi,epsilon) $turbulence;
div((nuEff*dev2(T(grad(U))))) Gauss linear;
}
laplacianSchemes
{
default Gauss linear corrected;
}
interpolationSchemes
{
default linear;
}
snGradSchemes
{
default corrected;
}
wallDist
{
method meshWave;
}
5.3 fvSolution文件
fvSolution文件指定线性方程组的求解方法。本算例中fvSolution文件如下所示。
FoamFile
{
version 2.0;
format ascii;
class dictionary;
object fvSolution;
}
// * * * * * * * * * * * * * * //
solvers
{
p
{
solver GAMG;
smoother GaussSeidel;
tolerance 1e-6;
relTol 0.1;
}
"(U|k|omega|epsilon)"
{
solver smoothSolver;
smoother symGaussSeidel;
tolerance 1e-6;
relTol 0.1;
}
}
SIMPLE
{
residualControl
{
p 1e-4;
U 1e-4;
"(k|omega|epsilon)" 1e-4;
}
nNonOrthogonalCorrectors 0;
pRefCell 0;
pRefValue 0;
}
potentialFlow
{
nNonOrthogonalCorrectors 10;
}
relaxationFactors
{
fields
{
p 0.3;
}
equations
{
U 0.7;
"(k|omega|epsilon).*" 0.7;
}
}
6 并行计算
前面演示的是采用串行计算,实际上对于稍微复杂的3D模型,建议使用并行计算。若想要采用并行计算,可以先往system
文件夹中添加字典文件decomposeParDict
,
foamGet decomposeParDict
并修改文件内容为:
FoamFile
{
version 2.0;
format ascii;
class dictionary;
object decomposeParDict;
}
// 设置使用40核进行计算
numberOfSubdomains 40;
// 设置并行方法为scotch,该方法无需任何参数
method scotch;
然后在终端中输入命令:
decomposePar
这样就可以将计算区域分解。然后计算可以使用下面的命令:
mpirun -np 40 simpleFoam -parallel
计算完毕后利用命令reconstruct
收集计算结果。
reconstructPar
然后就可以进行后处理了。
7 计算结果
-
压力分布
-
速度分布
-
离地10m横截面上速度分布
本篇文章来源于微信公众号: CFD之道
评论前必须登录!
注册