Pimple中对速度场的操作

ice_flow

各位老师好，

我在pimple的速度场中自己添加了一个扰动变量，编译可以通过，但是运行算例时会报错。
检查了很久不知道问题出在哪里，请大家帮忙。
这个是我在pimple算法的UEqn.H文件中添加的扰动UB：

// Solve the Momentum equation

MRF.correctBoundaryVelocity(U);

tmp<fvVectorMatrix> tUEqn
(
    fvm::ddt(U) + fvm::div(phiB, U)+ fvm::div(phi, UB)
  + MRF.DDt(U)
  + turbulence->divDevReff(U)
 ==
    fvOptions(U)
);
fvVectorMatrix& UEqn = tUEqn.ref();

单纯增加了一项变量，所以在主函数里面做了定义：

 // Define phiB and UB
    volVectorField UB
    (
        IOobject
        (
            "UB",
            runTime.timeName(),
            mesh,
            IOobject::READ_IF_PRESENT,
            IOobject::AUTO_WRITE
        ),
        mesh,
        dimensionedVector("UB", dimVelocity, vector::zero)
    );

    surfaceScalarField phiB
    (
        IOobject
        (
            "phiB",
            runTime.timeName(),
            mesh,
            IOobject::READ_IF_PRESENT,
            IOobject::AUTO_WRITE
        ),
        mesh
    );

编译是能够顺利通过的，运行算例的时候报错如下：

--> FOAM FATAL ERROR: (openfoam-2206)
Incompatible fields for operation
    [U] + [UB]

    From void Foam::checkMethod(const Foam::fvMatrix<Type>&, const Foam::fvMatrix<Type>&, const char*) [with Type = Foam::Vector<double>]
    in file /home/lb/v2206/OpenFOAM-v2206/src/finiteVolume/lnInclude/fvMatrix.C at line 1886.

FOAM aborting

#0  Foam::error::printStack(Foam::Ostream&) at ??:?
#1  Foam::error::simpleExit(int, bool) at ??:?
#2  void Foam::checkMethod<Foam::Vector<double> >(Foam::fvMatrix<Foam::Vector<double> > const&, Foam::fvMatrix<Foam::Vector<double> > const&, char const*) in ~/v2206/OpenFOAM-v2206/platforms/linux64GccDPInt32Opt/bin/dirpimpleFoam
#3  ? in ~/v2206/OpenFOAM-v2206/platforms/linux64GccDPInt32Opt/bin/dirpimpleFoam
#4  ? in ~/v2206/OpenFOAM-v2206/platforms/linux64GccDPInt32Opt/bin/dirpimpleFoam
#5  ? in /lib/x86_64-linux-gnu/libc.so.6
#6  __libc_start_main in /lib/x86_64-linux-gnu/libc.so.6
#7  ? in ~/v2206/OpenFOAM-v2206/platforms/linux64GccDPInt32Opt/bin/dirpimpleFoam
Aborted (core dumped)

提示场不匹配，所以我在循环的前后把两个场都打印了：

Starting time loop

UB field type: N4Foam14GeometricFieldINS_6VectorIdEENS_12fvPatchFieldENS_7volMeshEEE
U field type: N4Foam14GeometricFieldINS_6VectorIdEENS_12fvPatchFieldENS_7volMeshEEE
UB field dimension: [0 1 -1 0 0 0 0]
U field dimension: [0 1 -1 0 0 0 0]
Courant Number mean: 0 max: 0
Time = 0.005

UB field type: N4Foam14GeometricFieldINS_6VectorIdEENS_12fvPatchFieldENS_7volMeshEEE
U field type: N4Foam14GeometricFieldINS_6VectorIdEENS_12fvPatchFieldENS_7volMeshEEE
UB field dimension: [0 1 -1 0 0 0 0]
U field dimension: [0 1 -1 0 0 0 0]
PIMPLE: iteration 1

    From Foam::fvMatrix<Type>::fvMatrix(const Foam::GeometricField<Type, Foam::fvPatchField, Foam::volMesh>&, const Foam::dimensionSet&) [with Type = Foam::Vector<double>]
    in file /home/lb/v2206/OpenFOAM-v2206/src/finiteVolume/lnInclude/fvMatrix.C at line 360
    Constructing fvMatrix<Type> for field U

    From Foam::fvMatrix<Type>::fvMatrix(const Foam::GeometricField<Type, Foam::fvPatchField, Foam::volMesh>&, const Foam::dimensionSet&) [with Type = Foam::Vector<double>]
    in file /home/lb/v2206/OpenFOAM-v2206/src/finiteVolume/lnInclude/fvMatrix.C at line 360
    Constructing fvMatrix<Type> for field U

    From Foam::fvMatrix<Type>::fvMatrix(const Foam::GeometricField<Type, Foam::fvPatchField, Foam::volMesh>&, const Foam::dimensionSet&) [with Type = Foam::Vector<double>]
    in file /home/lb/v2206/OpenFOAM-v2206/src/finiteVolume/lnInclude/fvMatrix.C at line 360
    Constructing fvMatrix<Type> for field UB

    From Foam::fvMatrix<Type>::fvMatrix(const Foam::GeometricField<Type, Foam::fvPatchField, Foam::volMesh>&, const Foam::dimensionSet&) [with Type = Foam::Vector<double>]
    in file /home/lb/v2206/OpenFOAM-v2206/src/finiteVolume/lnInclude/fvMatrix.C at line 360
    Constructing fvMatrix<Type> for field U

    From Foam::fvMatrix<Type>::fvMatrix(const Foam::GeometricField<Type, Foam::fvPatchField, Foam::volMesh>&, const Foam::dimensionSet&) [with Type = Foam::Vector<double>]
    in file /home/lb/v2206/OpenFOAM-v2206/src/finiteVolume/lnInclude/fvMatrix.C at line 360
    Constructing fvMatrix<Type> for field U

    From Foam::fvMatrix<Type>::~fvMatrix() [with Type = Foam::Vector<double>]
    in file /home/lb/v2206/OpenFOAM-v2206/src/finiteVolume/lnInclude/fvMatrix.C at line 458
    Destroying fvMatrix<Type> for field U


--> FOAM FATAL ERROR: (openfoam-2206)
Incompatible fields for operation
    [U] + [UB]

    From void Foam::checkMethod(const Foam::fvMatrix<Type>&, const Foam::fvMatrix<Type>&, const char*) [with Type = Foam::Vector<double>]
    in file /home/lb/v2206/OpenFOAM-v2206/src/finiteVolume/lnInclude/fvMatrix.C at line 1886.

FOAM aborting

#0  Foam::error::printStack(Foam::Ostream&) at ??:?
#1  Foam::error::simpleExit(int, bool) at ??:?
#2  void Foam::checkMethod<Foam::Vector<double> >(Foam::fvMatrix<Foam::Vector<double> > const&, Foam::fvMatrix<Foam::Vector<double> > const&, char const*) in ~/v2206/OpenFOAM-v2206/platforms/linux64GccDPInt32Opt/bin/dirpimpleFoam
#3  ? in ~/v2206/OpenFOAM-v2206/platforms/linux64GccDPInt32Opt/bin/dirpimpleFoam
#4  ? in ~/v2206/OpenFOAM-v2206/platforms/linux64GccDPInt32Opt/bin/dirpimpleFoam
#5  ? in /lib/x86_64-linux-gnu/libc.so.6
#6  __libc_start_main in /lib/x86_64-linux-gnu/libc.so.6
#7  ? in ~/v2206/OpenFOAM-v2206/platforms/linux64GccDPInt32Opt/bin/dirpimpleFoam
Aborted (core dumped)

可以看到初始状态的两个场是类型一致的

UB field type: N4Foam14GeometricFieldINS_6VectorIdEENS_12fvPatchFieldENS_7volMeshEEE
U field type: N4Foam14GeometricFieldINS_6VectorIdEENS_12fvPatchFieldENS_7volMeshEEE
UB field dimension: [0 1 -1 0 0 0 0]
U field dimension: [0 1 -1 0 0 0 0]

但是进入第一个速度修正之后就报错了，说明这里出了问题。但是根据我自己的理解，后面的速度修正是通过添加网格运动，压力边界这些源项进行的，并没有改变速度的类型。或许变成矩阵了，但这些操作同时作用在两个场上，不知道问题出现在哪里。
后面我又用icofoam做了相同的修改，是可以运行的，感觉是在网格运动上吗，但是调试代码中也没有出现对应的报错，十分困扰。

李东岳

fvm::div(phi, UB)改成fvc::div(phi, UB)

ice_flow

可以了！谢谢老师。是因为这个fvm::div(phi, UB)迭代里只能有一个变量吗。

李东岳

你的方程是对U进行离散，fvm::div(phi, UB)还对UB进行离散。一个方程不能进行两个未知数的离散。

ice_flow

@李东岳明白了，谢谢老师！

Rachel0096

@ice_flow 您好，我也遇到了类似的问题，可以向您请教吗？这是对openfoam2.3.x版本中coalChemistryFoam的改写。大致也是在一个方程中，有两个不同的变量，编译可以通过，但在计算式出现报错incompatible fields for operation [N2] == [N2]。个人感觉也是同样的问题，但不是很清楚要如何去更改。YiEqn部分代码如下：

    for(label i=0; i<Y.size(); i++)
    {
        {
            // volScalarField Yi   = Ysf[nsf][i];
            volScalarField Ymean = Y[i];

            volScalarField tempYi
            (
                IOobject
                (
                    Y[i].name(),
                    runTime.timeName(),
                    mesh,
                    IOobject::NO_READ,
                    IOobject::NO_WRITE
                ),
                0*Y[0],
                "zeroGradient"
            );

            tempYi = Ysf[nsf][i];

            WienerTermYi = pow(2*gamma,0.5)*( fvc::grad(tempYi) & dW_sf);
            fvScalarMatrix YiEqn
            (
                  fvm::ddt(rho, tempYi)
                + mvConvection->fvmDiv(phi, tempYi)
                - fvm::laplacian(turbulence->muEff()/Sct, tempYi)
                ==
                  coalParcels.SYi(i, Ymean)
                + combustion->R(Ymean)
                + fvOptions(rho, Ymean)
		        + (rho*WienerTermYi/deltaT)
            );
            YiEqn.relax();
            fvOptions.constrain(YiEqn);
            YiEqn.solve(mesh.solver("tempYi"));
            YiEqn.solve();
            fvOptions.correct(tempYi);

            tempYi.max(0.0);
	        Ysf[nsf][i] = tempYi;
            Yt += tempYi;
        }
    }

李东岳

incompatible fields for operation [N2] == [N2]这是单位不一致。你植入的代码对应的方程单位不一致。你查查

李东岳

啊说错了，不是单位不一致，是赋值不统一。代码太多了，判断下哪一行有错误然后再进一步debug

ice_flow

@Rachel0096 这是多组分计算吗，我觉得你要检查下N2组分的初始定义，或者在程序里输出看一下是不是预想的值。我这里的问题是出现了多变量，应该不太一样。
类似这样

Info << "N2 field range: " << min(Y["N2"]).value() 
     << " to " << max(Y["N2"]).value() << endl;

Rachel0096

感谢老师们！这里确实是多变量的计算，类似于把单个组分场拆分为很多不同的场来计算，但是输入相同的源项，所以会有两种不同的Y在同一个方程里。程序就是在方程这一步出现了错误，我试试各位老师的方法，多输出一些内容看看哪里出问题了。

CFD中文网

Pimple中对速度场的操作