interFoam修改随时间变化的重力

danao

interFoam添加体积力，想通过修改重力，得到一个随时间和位置变化的重力。而在"readGravitationalAcceleration.H"里面，g的类型是 uniformDimensionedVectorField，而g与mesh.C()有关，所以编译时报错。
把类型改为DimensionedVectorField，报错：error: ‘DimensionedVectorField’ was not declared in this scope。改为volVectorField, 报错：‘Foam::tmp<Foam::Field<double> >’ is not derived from ‘const Foam::dimensioned<Type>’ ghf = (g & mesh.Cf()) - ghRef; ghRef这一项被标红了。应该改为什么类型才好？

Prometheus10

@danao 您好， danao！我也想设置一个随时间和位置变化的重力，您解决这个问题了吗？

李东岳

最简单的就是强制性的在不同的时间步下，改一下uniformDimensionedVectorField的g的值就行了？类似

while (piso.loop())
{
    g.value() = runTime.value(); 
}

Prometheus10

@李东岳（膜拜~）谢谢东岳老师！

Prometheus10

@李东岳
东岳老师您好！

因为对OpenFOAM的编程还不是很熟悉，目前我通过修改interIsoFoam做了三种尝试，想得到重力能随时间变化的效果，比如[g=sin(2 * M_PI * runTime.value)]。

过程中，我拷贝并修改~/application/solvers下的interIsoFoam文件夹后，将其命名为新的编译器，例如newinterIsoFoam。接下来陈述的三种方式都顺利编译成功，运行之前controlDoct我也有检查。

随后，我用damBreak案例来观察调整前后的情况，发现三种方式都没有使重力发生改变。依然是保持constant文件夹下的g [0, -9.81,0]，能否请您再出山，指点指点又菜又爱玩OpenFOAM的我，谢谢东岳老师!

我用的版本是dyfluid.com无私的OpenFOAM-v1812

第一种尝试：
按照您在本贴中的上一次提示，我将interIsoFoam.C的while (pimple.loop())下添加代码如下：

起始于interIsoFoam.C 110行~
        // --- Pressure-velocity PIMPLE corrector loop
                while (pimple.loop())
                {
            g.component(1).value() = (-1)*runTime.value(); //***新加内容***
                if (pimple.firstIter() || moveMeshOuterCorrectors)
.....

我从代码读出的效果应该是，Y轴负方向的值会随着pimple.loop的进行，g值从[0, 0, 0]变为[0, -2, 0]，但运行之后，依然是保持constant文件夹下的g [0, -9.81,0]。

第二种尝试：
因为v1812现在的readGravitationalAcceleration.H好像和数字版本 有一些不同

v版本的readGravitationalAcceleration.H是：

Info<< "\nReading g" << endl;
const meshObjects::gravity& g = meshObjects::gravity::New(runTime);

数字版本的readGravitationalAcceleration.H是：

     Info<< "\nReading g" << endl;
     uniformDimensionedVectorField g
     (
         IOobject
         (
             "g",
             runTime.constant(),
             mesh,
             IOobject::MUST_READ,
             IOobject::NO_WRITE
         )
     );

我将数字版本的readGravitationalAcceleration.H，复制到v1812上也能正常编译通过，通过您在另一篇帖子的回复--定义加速度的phi
我对照着做了以下修改：

     Info<< "\nReading g" << endl;
     uniformDimensionedVectorField g
     (
         IOobject
         (
             "g",
             runTime.constant(),
             mesh,
             IOobject::MUST_READ,
             IOobject::NO_WRITE
         );

    while (piso.loop())   //***新加内容***
        {   
             g.component(1).value() = (-1)*runTime.value(); //***新加内容***
        }  
    );

运行之后，damBreak依然是保持constant文件夹下的g [0, -9.81,0]

第三种尝试：
我照着另一份帖子说的按照改UEqn.H的方式引入附加体积力，改了interIsoFoam.C旁边的UEqn.H

    MRF.correctBoundaryVelocity(U);
    # define omega 2  //***新加内容***
    const dimensionedVector mySource("mySource", dimensionSet(1,-2,-2,0,0,0,0), 100*Foam::sin(M_PI*runTime.value()*omega)*vector(0,1,0)); //***新加内容***
//我粗略的考虑了水密度为1，因此直接把体积力的量纲用到了这里
    fvVectorMatrix UEqn
    (
        fvm::ddt(rho, U) + fvm::div(rhoPhi, U)
      + MRF.DDt(rho, U)
      + turbulence->divDevRhoReff(rho, U)
     ==
        fvOptions(rho, U)
      + mySource   //***新加内容***
    );

运行之后，damBreak依然是保持constant文件夹下的g [0, -9.81,0]。

实在是想了好久，三种方法都不行，当场心情就像溃坝一样damBreak!

请东岳老师出山！

xpqiu

@prometheus10
constant 下面的 g 文件，只是一个参数设置文件，默认情况下求解器会从这个文件中读取g 的值。你在程序里面修改了g 的值，除非你把修改后的值重新写出到 g 文件，否则 g 文件就不会变了。也就是说，你修改的重力是否生效了，应该从流场结果来判断，而不是通过g文件是否发生了变化来判断。

Prometheus10

@xpqiu xpqiux老师您好！
感谢您的回复，我现在才意识到了自己在这个环节的表述有点模糊。
我陈述的意思就如同您说的，。
我用上述三种方式改变重力之后，得到的流场结果依然还是g在 [0, -9.81,0]的流场结果，也就意味着我更改重力的三种方式并没有生效，所以一直很苦恼。

再次感谢您的回复！

CFD红皮书目前才看到第九章，编程还没开始看，所以一直不能找到问题在哪，希望您指点！

xpqiu

@prometheus10
哦，明白了。你的前两种方式，等于是都在尝试修改 g 这个变量的值，但是查看 interIsoFoam求解器的代码，它其实是共用了 interFoam 里面的 UEqn.H ，在这个文件里面，重力是通过如下这段代码起作用的，

            fvc::reconstruct
            (
                (
                    mixture.surfaceTensionForce()
                  - ghf*fvc::snGrad(rho)
                  - fvc::snGrad(p_rgh)
                ) * mesh.magSf()
            )

也就是说，真正进入到 UEqn 的是 ghf 这个变量。这个变量的定义在 gh.H 这个文件，如下：

    dimensionedScalar ghRef
    (
        mag(g.value()) > SMALL
      ? g & (cmptMag(g.value())/mag(g.value()))*hRef
      : dimensionedScalar("ghRef", g.dimensions()*dimLength, 0)
    );
    volScalarField gh("gh", (g & mesh.C()) - ghRef);
    surfaceScalarField ghf("ghf", (g & mesh.Cf()) - ghRef);

所以这就是为啥你修改 g 的值不起作用了。因为g 的值只是在开始阶段，读取进来，并据此计算了 ghf 的值，然后 ghf 进入 UEqn。要想让重力发生变化，需要 ghf进入 UEqn之前修改其值才行。或者，在计算 ghf 之前，修改 g 的值。

至于第三种体积力的方式为啥不起作用，我没看出来，主要是不清楚你说的改了 "改了interIsoFoam.C旁边的UEqn.H" 具体指的是什么，你是直接改了 interFoam 那边的 UEqn.H，然后重新编译的 interIsoFoam？

Prometheus10

@xpqiu xpqiu老师您好！感谢您的回复，我大致明白了您对前两个方式的解释！谢谢您！

至于在第三种方式中，提到的 "改了interIsoFoam.C旁边的UEqn.H"
是指的UEqn.H和interIsoFoam.C文件都在同一目录下，即：~/applications/solvers/multiphase/interIsoFoam目录。

以后发帖我一定多多注意自己的措辞，耽误了您的时间，同时为上一次打错了您的昵称抱歉。
再次谢谢您xpqiu老师