《icoFoam解析》中压力修正的代码问题

[jzheng]

请教各位大神个问题，在《icoFoam解析》中求第一次修正压力场pr的公式是~~~~

但是，在pEqn的代码是

我的问题是，把fvc::div(phiHbyA)代码还原成公式的话是∇∙(phiHbyA)，它是怎么和修正公式里面的∇∙(HbyA)对应上的呢？

[李东岳]

@jzheng 说：

∇∙(phiHbyA)

Hi，

HbyA是矢量，phiHbyA是面标量。fvc::div(phiHbyA)表示的就是$\nabla \cdot HbyA$。再比如fvc::div(phi)就是$\nabla \cdot U$。这个div被操作符重载了。

[散漫守望2016]

竟然和我有一样的疑惑，当时这里一直不理解，现在还是不太理解。东岳兄那个乱码了，可以改改。

[jzheng]

@李东岳 明白了，谢谢东岳兄，厉害

[Wayne]

补充一下。

公式乱码刷新一下就行。

fvc::div代码如下：

template<class Type>
tmp<GeometricField<Type, fvPatchField, volMesh> >
div
(
    const GeometricField<Type, fvsPatchField, surfaceMesh>& ssf
)
{
    return tmp<GeometricField<Type, fvPatchField, volMesh> >
    (
        new GeometricField<Type, fvPatchField, volMesh>
        (
            "div("+ssf.name()+')',
            fvc::surfaceIntegrate(ssf)
        )
    );
}

fvc::div(phiHbyA) 就是对单元体表面的phiHbyA进行求和。这个依据高斯定理，相当于对 ∇∙(HbyA) 在单元体上进行积分。

所以这个问题的关键是，div和laplacian等操作符是对
\begin{equation}
\nabla \cdot (HbyA^r) = \nabla \cdot(\frac{1}{A_{\mathrm{p},f}} \nabla p^r)
\end{equation}
这个半离散方程进行进一步的离散，即积分。

[Wayne]

另外，
参考
fvc::div(u)=fvc::div(phi)。

对应的代码如下：

template<class Type>
tmp
<
    GeometricField
    <typename innerProduct<vector, Type>::type, fvPatchField, volMesh>
>
gaussDivScheme<Type>::fvcDiv
(
    const GeometricField<Type, fvPatchField, volMesh>& vf
)
{
    tmp
    <
        GeometricField
        <typename innerProduct<vector, Type>::type, fvPatchField, volMesh>
    > tDiv
    (
        fvc::surfaceIntegrate
        (
            this->mesh_.Sf() & this->tinterpScheme_().interpolate(vf)
        )
    );

    tDiv().rename("div(" + vf.name() + ')');

    return tDiv;
}

即div(U)的计算方法也是先插值到单元体表面上，然后进行求和。

[李东岳]

@Wayne

Updated. 非常感谢。

[张某人]

我找到了一个对div和snGrad的理解，希望可以帮上忙，内容来源：Giskard's CFD Learning Tricks
https://wenku.baidu.com/view/838580050975f46527d3e1ff.html
我的理解是fvc::div(phiHbyA)已经将原控制方程中的散度表示为高斯积分的求和形式了，不知道这样理解是否正确？

[李东岳]

我的理解是fvc::div(phiHbyA)已经将原控制方程中的散度表示为高斯积分的求和形式了

是的。正确。