OpenFOAM动网格的壁面边界通量修正

OpenFOAM使用动网格时，壁面速度边界条件要换成movingWallVelocity类型。这个边界类型是如何实现壁面运动时通量的计算和修正的？上一篇学习了内场网格变化时通量的修正，这一篇学习这一movingWallVelocity类型的实现原理。

movingWallVelocity边界条件类型继承自fixedValue边界类型。文件在src/finiteVolume/fields/fvPatchFields/movingWallVelocity中找到。此类型是可以向下兼容的。如果网格没有变形，那么会退化到fixedValue边界类型。实现代码如下：

void Foam::movingWallVelocityFvPatchVectorField::updateCoeffs()
{
    if (updated())
    {
        return;
    }

    const fvMesh& mesh = internalField().mesh();

    if (mesh.moving())
    {
        const fvPatch& p = patch();
        const polyPatch& pp = p.patch();
        const pointField& oldPoints = mesh.oldPoints();

        vectorField oldFc(pp.size());

        forAll(oldFc, i)
        {
            oldFc[i] = pp[i].centre(oldPoints);
        }

        const scalar deltaT = mesh.time().deltaTValue();

        const vectorField Up((pp.faceCentres() - oldFc)/deltaT);

        const volVectorField& U =
            static_cast<const volVectorField&>(internalField());

        scalarField phip
        (
            p.patchField<surfaceScalarField, scalar>(fvc::meshPhi(U))
        );

        const vectorField n(p.nf());
        const scalarField& magSf = p.magSf();
        tmp<scalarField> Un = phip/(magSf + VSMALL);


        vectorField::operator=(Up + n*(Un - (n & Up)));
    }

    fixedValueFvPatchVectorField::updateCoeffs();
}

速度在壁面上的边界值，通过movingWallVelocity计算得到，在U.correctBoundaryConditions()调用。可以看到Up就是一开始patch的速度值，用patch在网格变形前后的变化率计算得到（const vectorField Up((pp.faceCentres() - oldFc)/deltaT);）。然后，通过meshPhi计算到phip，即网格关于U的通量，与时间格式有关。精度最多到二阶，能提高网格运动时通量修正的精度。Un相当于interpolate(U)，即面上速度。也就是说，Up的计算涉及网格空间的变化，那么之后还要在时间变化上再修正一次。可以看到，最后输出（Up + n*(Un - (n & Up)）。第二项可以认为是时间精度格式带来的修正。关于meshPhi可以看我上一篇推文：OpenFOAM动网格的通量修正_陈十七的博客-CSDN博客