CFD中文网

@yangom 无滑移边界条件，下面就是固壁面，无滑移就是壁面上速度为0，那下表面就需要从入口处开始发展发展出一个边界层。在这个边界层中，沿壁面法向向外，速度当然是从0逐渐变化成外流速度的
你直接看也看不出来有没有速度，画几条曲线出来。
1 靠近近壁面某个高度上，沿流动方向的一条线，画这条线上沿流动方向的速度
2 选几个位置，沿垂直于壁面方向的一条线，画这条线上沿流动方向的速度

小球这里没速度是和你速度边界条件的设置有关系，你应该是把下面那个边界设置成了noSlip。
需要更多的信息。

加源项的时候，这个源项的表达式好像有点问题，比如下面添加vectorSemiImplicitSource的时候源项是( (299.4600 0 0) 0)而不是(299.4600 0 0)。原因是加入的源项表达式为：
$$
S=S_C+S_Pϕ_P
$$
标量源项的代码也是这么写的：

        injectionRateSuSp
        {
            variable_name     (Sc Sp);
        }

momentumSource
{
   type vectorSemiImplicitSource;
   active on;
   

   vectorSemiImplicitSourceCoeffs
   {
      selectionMode all;
      //volumeMode      absolute; // specific
      volumeMode        specific;
      injectionRateSuSp
      {
         U           ( (299.4600 0 0) 0); //partial p / partial x
      }
   }
}

两个参考：
https://caefn.com/openfoam/fvoptions-acousticdampingsource
https://xiaopingqiu.github.io/2016/03/20/fvOptions2/

刚又找了一下之前的笔记，不是上面这个原因，因为你用的是vectorCodedSource。我之前也刚好写过scalarCodedSource，也发现量纲不对导致的问题。之前记录的笔记是：
有一个问题是感觉并不需要乘以cell的体积。官方给的文档里面是有乘以体积的，但这个网页就很迷惑，https://www.openfoam.com/documentation/guides/v2012/doc/guide-fvoptions-sources-coded.html ，上面的表达式是除以体积，下面的代码又是乘以体积了。找到一个人的注释 https://xiaopingqiu.github.io/2016/03/20/fvOptions2/

// fvMatrix<Type> 类中对“+=”操作符进行了重载，所以，eqn与Su的相加，相当于eqn+Su*mesh.V()，要不然eqn与Su的量纲不一致。
eqn += Su + fvm::SuSp(Sp, psi);

所以还是要乘以体积。

有一个在区域当中加源项以实现消波的模块
https://caefn.com/openfoam/fvoptions-acousticdampingsource

function object surfaces
实时输出就用

writeControl timeStep;
writeInterval 1;

@SSSSK 这两个图横纵坐标分别是什么？看起来确实是从20s开始统计的话收敛很快，150s开始统计收敛就好慢了，而且150-250峰值也没到2.5以上。有参考的DNS数据吗？和DNS比较一下。
你是用的哪一种驱动槽道的方式？是meanVelocityForce的话，就把pressureGradient随时间变化的图像贴一下。

你画一下
20s开始统计，t=21、22……30s 的UPrime2Mean；
和
100s开始统计，t=101、102……110s 的UPrime2Mean。

@SSSSK
restartTime 20;
和
startTime 20;
在20s之后输出的UMean，是一样的。
湍流是否完全发展，还可以看摩擦雷诺数的时间历程；10、20、30……这样不同时刻的速度剖面。得综合判断。

我感觉是不是从100s开始统计的时候，其实是把前面的时间步也考虑了，在平均下来就会很小

100s开始统计不会考虑前面时间步的。你想如果它能一瞬间把前面100s的数据都考虑进去的话，那就是前面100s内每个时间步的全流场数据都存着，没那么大的内存或者硬盘容量的。

研究一下源代码吧。有个大佬的笔记讲得挺好↓
https://xiaopingqiu.github.io/2015/04/12/fieldAverage/

我用的不深入。OpenFOAM是有这种边界条件的：

bottomWall
    {
        type            constantAlphaContactAngle;
        theta0          0;
        limit           gradient;
        value           uniform 0;
    }

我也一直纠结这玩意儿。

如果你算channel的时候用的是meanVelocityForce，那么程序会根据UMean来调整压力梯度。当你在20s开始统计的时候，最初的几个时间步UMean波动会很大，连带的添加的这个压力梯度也会很大幅度地变化，流场还需要很久才能平稳。

momentumSource
{
    type            meanVelocityForce;
    active          yes;
 
    meanVelocityForceCoeffs
    {
        selectionMode   all;
 
        fields          (U);
        Ubar            (1.0363 0 0);
        relaxation      1.0;
    }
}

就算不用meanVelocityForce，用的是vectorSemiImplicitSource。在重启FieldAverage的时候，最初的一段时间UMean还是波动很大的，UPrime2Mean也需要一段时间才能平稳。

最根本的原因还是FieldAverage它每一步计算平均值，都是将当前时间步的U和上一个时间步存下来的UMean进行计算，只有算一段时间后，才能平稳。

下面是我用的v2012版本的FieldAverage介绍，虽然这里是几个数的平均值，但实际上计算的时候也会将时间、时间步长考虑在内。但要知道的就是这个平均值每一个时间步都是新算出来的，平方均值也是每一个时间步新算出来的。

你20s发展完全，使用FieldAverage统计可能还得往后面算。或者自己以一定时间间隔保存全场数据然后进行后处理。

之前学习的时候看到过这个文章，挺详细的，分享一下
https://blog.csdn.net/hanbingchegu/article/details/107160364

最近在写论文，因为用的是interFoam求解器，所以在写控制方程的时候一致参考李东岳老师的这个笔记 http://dyfluid.com/interFoam.html 。最近仔细看了下密度、运动粘度、动力粘度之间的关系，发现一直以来都没注意到的错误：
错误在于公式24

$$
\nabla\cdot\boldsymbol{\tau}=\nabla\cdot\left(\nu\left(\nabla\bfU+\nabla\bfU^{\mathbf{T}}\right)\right)
$$

在多相流的控制方程里面，密度$\rho$是和速度放一起了，所以上面公式当中的 $\nu$ 应该为 $\mu$ 。

interFOAM类的求解器用的是不可压缩不相容的混合物模型，在代码中是这样定义的

immiscibleIncompressibleTwoPhaseMixture mixture(U, phi);

它要求输入每一个相的密度$\rho$和运动粘度$\nu$，而$\nu$是在$\mu$的基础上计算的。
在immiscibleIncompressibleTwoPhaseMixture 这个类的父类incompressibleTwoPhaseMixture中有计算类成员nu_的函数：

void Foam::incompressibleTwoPhaseMixture::calcNu()
{
    nuModel1_->correct();
    nuModel2_->correct();

    const volScalarField limitedAlpha1
    (
        "limitedAlpha1",
        min(max(alpha1_, scalar(0)), scalar(1))
    );

    // Average kinematic viscosity calculated from dynamic viscosity
    nu_ = mu()/(limitedAlpha1*rho1_ + (scalar(1) - limitedAlpha1)*rho2_);
}

而$\mu$的计算公式为 $\mu = \alpha_1 \rho_1 \nu_1 + (1-\alpha_1) \rho_2 \nu_2$：

Foam::tmp<Foam::volScalarField>
Foam::incompressibleTwoPhaseMixture::mu() const
{
    const volScalarField limitedAlpha1
    (
        min(max(alpha1_, scalar(0)), scalar(1))
    );

    return tmp<volScalarField>::New
    (
        "mu",
        limitedAlpha1*rho1_*nuModel1_->nu()
      + (scalar(1) - limitedAlpha1)*rho2_*nuModel2_->nu()
    );
}

也就是说$\nu$公式和东岳老师笔记 http://dyfluid.com/interFoam.html 的公式25不一致，正确的公式是
$$\nu = \frac{{{\alpha _1}{\rho _1}{\nu _1} + \left( {1 - {\alpha _1}} \right){\rho _2}{\nu _2}}}{{{\alpha _1}{\rho _1} + \left( {1 - {\alpha _1}} \right){\rho _2}}}$$

那么粘性应力在程序中是怎样计算的？
interFOAM的动量方程在UEqn.H当中，是

fvVectorMatrix UEqn
    (
        fvm::ddt(rho, U) + fvm::div(rhoPhi, U)
      + MRF.DDt(rho, U)
      + turbulence->divDevRhoReff(rho, U)
     ==
        fvOptions(rho, U)
    );

可以看到粘性应力项和湍流粘度项混在一起，直接由湍流模型的对象turbulence给出。但是只有当给出具体湍流模型的时候，才能知道究竟是怎样的表达式。
我再往下去就有点晕了，不知道要找哪个源代码文件当中的具体函数了，求大佬指点一下。就是没找到

不去探求源代码的话，我的理解是，湍流模型给出的都是nut，在各个算例文件当中都给的是动力湍流粘度。
来自湍流的动力粘度为$\rho \nu_t$是要用当地的密度乘以湍流模型给出的运动粘度。
也就是说：这里的turbulence->divDevRhoReff(rho, U)在WALE模型下，或者其它的一些模型下，计算 粘性应力 + 湍流应力 在方程中的那一项，对应的表达式为
$$\nabla\cdot\boldsymbol{\tau}=\nabla\cdot\left( \mu_{eff} \left(\nabla\bfU+\nabla\bfU^{\mathbf{T}}\right)\right)$$
其中
$${\mu _{eff}} = {\mu} + {\mu _t} = {\mu} + \rho {\nu _t} = {\alpha _1}{\rho _1}{\nu _1} + \left( {1 - {\alpha _1}} \right){\rho _2}{\nu _2} + \left[ {{\alpha _1}{\rho _1} + \left( {1 - {\alpha _1}} \right){\rho _2}} \right]{\nu _t}$$

@李东岳 谢谢李老师关心呀。算例模拟完了，文章😂我还在写。
算出来的结果不是很好，气泡槽道流收敛性还比较差，几个参数之间的规律不是很明显，故事不太好讲。

可能计算的时间不够长，继续计算一下

frontAndBack empty这种边界条件是针对二维网格的吧，三维的网格应该不能用。
你试试把U和p的frontAndBack都改成zeroGradient

压力边界条件也要合适：出口压力设置成固定值，两个入口的压力设置成zeroGradient

就我粗浅的经验来看，网格质量太差。网格尺寸变化的地方，就是粗细网格交界的那里，粗网格比细网格大太多了，最好是渐变一些。

Shihang-Chen 在 cfdem中cfd流场重组reconstructPar 中说：

如果开自适应网格，就先
reconstructParMesh
然后
reconstructPar -noLagrangian
如果没开自适应
就直接reconstructPar -noLagrangian

看看这个。你是开自适应网格了吗？

else后面不要跟括号的表达式，去看一下C++的格式哇