[simpleFoam]、雷诺应力及湍流模型原理

[evensun]

各位老师好，我是一个刚接触湍流模型的萌新，关于simpleFoam和k-epsilon湍流模型有以下问题想请教：
①k-epsilon模型是通过k和epsilon两方程解出涡粘系数nut，进而求解雷诺应力使方程封闭对吗，雷诺应力是包含在turbulence->divDevReff(U)这项当中的吗；
②在simpleFoam中，UEqn中有turbulence->divDevReff(U)这一项，为什么我在src/TurbulenceModels文件夹下搜索“divDevReff”关键词却没有找到其表达式，只有divDevRhoReff的表达式呢；
③以下是我找到的nuEff的表达式，其有nu和nut组成，这是否是湍流模型在求解器计算UEqn时发挥作用的体现呢（不然SIMPLEC算法直接就解出U和p了）

//- Return the effective viscosity
        virtual tmp<volScalarField> nuEff() const
        {
            return volScalarField::New
            (
                IOobject::groupName("nuEff", this->alphaRhoPhi_.group()),
                this->nut() + this->nu()
            );
        }

④如果在求解器中设置了很多速度，比如U、Ut、Uf等，那么湍流模型在计算的时候，是选择其中的U来计算雷诺应力吗，选择的依据是什么呢；
我的湍流模型基础较薄弱，希望老师多多指教，感激不尽

[evensun]

谢谢各位老师，我找到了相关资料：
对于线性涡粘模型，divDevReff(U) = divDevRhoReff(U)
这样问题②就明白啦

[evensun]

不好意思各位老师，以上问题我已经全部明白了，以后我会多读代码加深理解的

[羽之下]

@evensun

[cccrrryyy]

@evensun 建议搞明白自己的问题后稍微辛苦下把自己的回答也写下来，这样以后别人看到了会有所帮助

[evensun]

雷诺应力确实有一部分包含在turbulence->divDevReff(U)这一项中，在OpenFoam中这一项称为有效雷诺应力，表达式如下
$$\begin{array}{}\text{(1)}& {\mathbf{R}}_{\mathrm{e}\mathrm{f}\mathrm{f}}={\nu }_{\mathrm{e}\mathrm{f}\mathrm{f}}[\mathrm{\nabla }\overline{\mathbf{u}}+(\mathrm{\nabla }\overline{\mathbf{u}}{)}^{\mathsf{T}}-\frac{2}{3}(\mathrm{\nabla }\cdot \overline{\mathbf{u}})\mathbf{I}]\end{array}$$其中的${\nu }_{\mathrm{e}\mathrm{f}\mathrm{f}}=\nu +{\nu }_t$，${\nu }_t$即为通过k-epsilon模型计算出的涡粘性系数。turbulence->divDevReff(U)在动量方程中，由此湍流模型在求解速度时发挥作用。
另外关于turbulence->divDevReff(U)这一项的调用，虽然此项在不可压求解器中，但如果采用线性涡粘的湍流模型，调用的是linearViscousStress 中的 divDevRhoReff(U)（可压缩也用这个），具体代码为：

template<class BasicTurbulenceModel>
Foam::tmp<Foam::fvVectorMatrix>
Foam::linearViscousStress<BasicTurbulenceModel>::divDevRhoReff
(
    volVectorField& U
) const
{
    return
    (
      - fvc::div((this->alpha_*this->rho_*this->nuEff())*dev2(T(fvc::grad(U))))
      - fvm::laplacian(this->alpha_*this->rho_*this->nuEff(), U)
    );
}

计算时根据括号中的量选择带入计算的是哪个速度。（其中负号是因为定义在系数矩阵等式的右边）
以上内容整理自OpenFOAM 中的不可压和可压 RANS 方程，网页内容非常详细，要是我早点看到能省不少麻烦

[evensun]

@cccrrryyy 谢谢老师提醒，我把我理解的内容整理在楼下了

[cccrrryyy]

感谢，辛苦啦！

[李东岳]

@evensun 感谢分享。

你用的是哪个雷诺应力模型？LRR？SSG？

[evensun]

@李东岳 我用的是RAS下的k-epsilon。其实目前的工作是在openfoam中植入新的湍流模型，但想通过学习k-epsilon的求解过程来打开思路。

[AppleKiller]

@evensun 大佬您好，这段代码中的this->alpha_表示什么呢？因为我感觉
this->rho_*this->nuEff())*dev2(T(fvc::grad(U))
已经足够表示
$${\mu }^{eff}(\mathrm{\nabla }UU+\mathrm{\nabla }{U}^{\mathrm{\top }}{U}^{\mathrm{\top }}-\frac{2}{3}\mathrm{\nabla }\cdot UU\mathcal{I})$$

[AppleKiller]

问题解决，this->alpha_的定义在

template
<
    class Alpha,
    class Rho,
    class BasicTurbulenceModel,
    class TransportModel
>
class TurbulenceModel
:
    public BasicTurbulenceModel
{
 
public:
 
    typedef Alpha alphaField;
    typedef Rho rhoField;
    typedef TransportModel transportModel;
 
 
protected:
 
    // Protected data
 
        const alphaField& alpha_;
        const transportModel& transport_;
     ..........
}

这里的alpha_应该表示相分数。