代码实现？

[小龙]

@李东岳 东岳老师

volVectorField unknown("unknown", fvc::div(tauU));

这个是不是应该改成

volScalarField unknown("unknown", fvc::div(tauU));

这一项应该是个标量场？

[李东岳]

tauU是二阶张量，做div后是矢量

[小龙]

@李东岳 这里的tau是μgrad(U)吗？（μ是动力粘度），还是可以在OpenFoam里调用的？

[李东岳]

tau是$\tau$，需要用Bousinessq假定来计算，参考 http://47.93.231.95/doku.php#boussinesq_approximation

[Izumi]

@李东岳 请问tau如下定义可以吗？

volTensorField tau("tau", muEff*(dev2(gradU.T()) + gradU));

因为
$$\begin{array}{}\text{(1)}& dev2((\mathrm{\nabla }u{)}^T)=(\mathrm{\nabla }u{)}^T-\frac{2}{3}tr((\mathrm{\nabla }u{)}^T)I=(\mathrm{\nabla }u{)}^T-\frac{2}{3}(\mathrm{\nabla }.u)I\end{array}$$
所以
$$\begin{array}{}\text{(2)}& \tau =\mu (dev2((\mathrm{\nabla }u{)}^T)+\mathrm{\nabla }u)\end{array}$$

[小龙]

@izumi tau应该是这样定义的

[程迪]

@izumi
不建议这么实现，参考rhoCentralFoam的实现，应该是把转置项显式离散，非转置项按laplacian隐式离散。

因为OF本身的限制，转置项会有各个分量间的耦合，所以必须显式离散，非转置项不用。

隐式的部分越多，时间步限制通常也越弱一些，所以OF一般实践上喜欢把粘性项拆拆拆之后尽可能多地采用隐式离散。实在不能隐式离散的源项、交叉项才用显式离散。

参考：https://chengdi123000.github.io/2018/01/05/rhoCentralFoam解析/#粘性应力项

[Izumi]

@程迪
谢谢！看rhoCentralFoam时没看懂，看了你的链接清楚一些了。
我主要是想在能量方程中添加$\mathrm{\nabla }\cdot (\tau \cdot u)$，然后计算耗散函数$\tau :(\frac{1}{2}(\mathrm{\nabla }u+(\mathrm{\nabla }u{)}^T))$。
根据CFD中的能量方程，$\mathrm{\nabla }\cdot (\tau \cdot u)$好像没法展开成可以隐式离散的项。
看rhoCentralFoam中的能量方程，除时间项外，用的都是显式离散。所以，对于能量方程来说，是不是不必拆项并使用隐式格式？
至于动量方程中的$\mathrm{\nabla }\cdot \tau$，对于一般的求解器，应该是包含在UEqn.H里的turbulence->divDevRhoReff(U)吧，不用我们自己定义。

rhoCentralFoam中的能量方程：

surfaceScalarField sigmaDotU
        (
            "sigmaDotU",
            (
                fvc::interpolate(muEff)*mesh.magSf()*fvc::snGrad(U)
              + fvc::dotInterpolate(mesh.Sf(), tauMC)
            )
          & (a_pos*U_pos + a_neg*U_neg)
        );

        solve
        (
            fvm::ddt(rhoE)
          + fvc::div(phiEp)
          - fvc::div(sigmaDotU)
        );

[程迪]

@izumi
我的理解是，拆项的话你得把phiEp和sigmaDotU表示成rhoE的函数才行，这似乎有点儿困难吧。。

[Izumi]

@程迪
rhoE应该是$\frac{\mathrm{\partial }\rho e}{\mathrm{\partial }t}+\frac{\mathrm{\partial }\rho K}{\mathrm{\partial }t}$
phiEp应该是$\mathrm{\nabla }\cdot (\rho \mathbf{U}e)+\mathrm{\nabla }\cdot (\rho \mathbf{U}K)+\mathrm{\nabla }\cdot (p\mathbf{U})$
sigmaDotU应该是$\mathrm{\nabla }\cdot (\tau \cdot \mathbf{U})$
感觉拆不了。

另外，请问对定常的求解器也有时间步限制吗？如果没有的话，是不是定常的求解器中，比如buoyanySimpleFoam，就不用管隐式还是显式了？