icoFoam UxIter_0=1000

Haining LUO

目标算例：基于圆管直径的雷诺数2400 (0.3 m/s * 8e-3 m/ 1e-6 m^/s)，没有湍流模型(icoFoam)
网格：结构网格，第一层网格yPlus~1（基于Re=2400边界层近似公式得到uTau [~0.02m/s]，yPlus=uTau*deltaY/nu），圆管中心的yPlus～80，径向网格数～30（是不是太稀疏？）
边界条件：
速度：两个进口（垂直）均为turbulentInlet，0.3m/s的均匀速度加上随机脉动，出口法向零梯度，壁面无滑移
压力：法向梯度为零，出口为压力uniform参照值

主要遇到两个问题：

很贴近边界[第二第三个]上的网格有较大的速度，与进口速度一个量级；靠近mixing zone能看到较强的速度场不连续性
考虑到网格可能太稀疏，有加倍网格（三个方向都几乎加倍），但无论怎么调小时间步长Ux Uy Uz 迭代(UxIter_0)第一个时间步都飙到了1000，最终残差下不了1e-2

对于网格质量(pointwise做的)，这是较密的网格，fvScheme 里面 nonOrthorgonal Corrector =1：

Time = 0

Mesh stats
    points:           779916
    faces:            2299824
    internal faces:   2262096
    cells:            760320
    faces per cell:   6
    boundary patches: 4
    point zones:      0
    face zones:       0
    cell zones:       0

Overall number of cells of each type:
    hexahedra:     760320
    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                  
    Port1               1056     1081     ok (non-closed singly connected)  
    Port2               1056     1081     ok (non-closed singly connected)  
    Port3               1056     1081     ok (non-closed singly connected)  
    wall                34560    34631    ok (non-closed singly connected)  

Checking geometry...
    Overall domain bounding box (-0.08 -0.08 -0.004) (0.08 0.004 0.004)
    Mesh (non-empty, non-wedge) directions (1 1 1)
    Mesh (non-empty) directions (1 1 1)
    Boundary openness (8.46492e-16 1.14034e-16 -8.71928e-17) OK.
    Max cell openness = 8.80197e-16 OK.
    **Max aspect ratio = 9.65639 OK.**
    Minimum face area = 8.44179e-09. Maximum face area = 4.11244e-07.  Face area magnitudes OK.
    Min volume = 2.00787e-12. Max volume = 4.4809e-11.  Total volume = 1.18509e-05.  Cell volumes OK.
    **Mesh non-orthogonality Max: 55.6104 average: 19.4413**
    Non-orthogonality check OK.
    Face pyramids OK.
    **Max skewness = 2.05024** OK.
    Coupled point location match (average 0) OK.

Mesh OK.

End

在较疏的case，我把nu变大十倍，nu=1e-5，UxIters_0 = {1,2,3}，最终残差可以收敛到1e-5，如下图（UMean在10*(D/U)时间的截面云图），进口是左边和下面，转角处有一个速度较大的网格点（离壁面第二个网格）。我觉得收敛得快是因为粘性项大了，数值稳定性强，所以我应该在雷诺数2400的算例改改fvScheme?还是得用其他的smooth solver?个人猜想速度大的那个网格产上的可能原因是我分块的时候网格随着分的块形变（下面看得很清楚），所以会有横向的flux，横向的flux通过什么输运到下游呢？通过那个斜着45度角的flux，但是仔细一看那个地方速度大的方向不是X而是Y方向（图中只有magnitude）。那会是什么原因呢？

附上实验数据，在那个角落速度不会那么大

所以问题是：

局部速度过大是个数值上的问题吗还是网格策略不对？
大雷诺数迭代不收敛应该是数值问题，是不是该试一试迎风格式什么的？还是换个U的求解器？

附上 fvScheme

// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //

ddtSchemes
{
    default         Euler;
}

gradSchemes
{
    default         Gauss linear;
    grad(p)         Gauss linear;
}

divSchemes
{
    default         none;
    div(phi,U)      Gauss linear;
}

laplacianSchemes
{
    default         Gauss linear orthogonal;
}

interpolationSchemes
{
    default         linear;
}

snGradSchemes
{
    default         orthogonal;
}

fluxRequired
{
    default         no;
    p               ;
}

李东岳

是说这个地方不符合实验么？

Haining LUO

@haining-luo
忘了加上：圆管截面上是O-H型的网格，能够贴体加密，算圆管流动没有问题

还有fvSolution

solvers
{
    p
    {
        solver          PCG;
        preconditioner  DIC;
        tolerance       1e-06;
        relTol          0;
    }

    U
    {
        solver          smoothSolver;
        smoother        symGaussSeidel;
        tolerance       1e-05;
        relTol          0;
    }
}

PISO
{
    nCorrectors     2;
    nNonOrthogonalCorrectors 1;
    pRefCell        0;
    pRefValue       0;
}

另外还有一个问题，因为问题是强三维的，所以N-S方程icoFoam segerated的方式（即依次求解Ux, Uy, Uz）的方式是不是不够好？(https://openfoam.org/release/2-2-0/matrix-solvers/)

我这个例子里面是否可以考虑改求解器顺序(Uy, Ux, Uz)来验证问题是不是强三维coupled，但应该怎么改呢？

Haining LUO

@李东岳对的，实验上无量纲化了，我算的无量纲后的速度是4，实验在那个位置小于2，贴近壁面的地方<0.25

李东岳

你贴个那个paraview的图看起来有点震荡，你这个网格多大？实验那个用PIV测的？icoFoam是层流求解器，但是你用的是湍流进口？

你可以试试他给的设置：

U
    {
        type            coupled;  // optional, defaults to segregated
        solver          PBiCCCG;
        preconditioner  DILU;
        tolerance       (1e-05 1e-05 1e-05);
        relTol          (0 0 0);
    }

不过我没试过这个。网格数量不大的话，不保密的话发给我看看我算算？你这个挺有意思，因为有实验验证。

Haining LUO

@李东岳

我用了上面的设置，正在算，同样的网格（Hexahedra: 23万，对应上面“稀疏”的case）coupled的solver收敛迅速！虽然如你所说震荡仍然很大，但迭代次数U[x,y,z]Iters_0=1000 => U[x,y,z]Iters_0={1,2,3}。看来对于我的算例symGaussSeidel明显不如PBiCCCG，一查结果发现后者解非对称而symGaussSeidel是解对称矩阵的。请问有办法知道自己求解的问题是否是对称的吗，比如检查A矩阵UEqn.A()？

图对应的网格很稀疏（Hexahedra: 23万），实验是PIV，是中间截面的数据

导师想让我做DNS，icoFoam应该可以解DNS，不过不知道PISO的DNS能力相比pimple怎么样？

我把Re2400_mesh2_coupledSolver放在网盘了：（这是稀疏的网格，用pointwise后来发现这个网格有有些手误的地方，应该做镜像然后merge，这样就完全对称；你可以也试一试，两个solver就是两重天，我这个白痴:crying: ）
https://pan.baidu.com/s/1smtVW1f
密码：w39p

李东岳

我正在计算还没结果。

不过从你的网格来看，能不能提高些网格质量，你给我的这个算例网格不是很完美。另外，由于你这几个几何结构比较简单，你可以试试snappyHexMesh做全正交网格。如果你用icoFoam，算法上很简单了并且没有什么模型，你只能提高网格质量来提高计算结果。

做DNS我建议用PISO，PIMPLE有松弛。

这些地方网格非正交很差，你还比较关心这个地方，但是也不太好处理，只能用snappyHexMesh搞一搞

李东岳

计算了几个时间步，我感觉你应该看看提高网格质量。如果你不关心这个圆圈部位的数据，应该也没问题。

CFD中文网

icoFoam UxIter_0=1000