1. Home
  2. OpenFOAM
  3. SprayFoam 只喷固体该如何设置？

SprayFoam 只喷固体该如何设置？

  • chpjz0391

    0_1458720018133_New Microsoft PowerPoint Presentation.jpg
    模型如图所示，我需要模拟一个辆车的尾气颗粒在一个隧道里怎么分布，颗粒与颗粒的碰撞需要考虑，这辆车只喷颗粒物质一定速度喷（比如每秒20000个），然后隧道一段以一定的速度进风。 之前@散漫守望2016 给我介绍可以看看coalchemistryFoam，但是我今天看了sprayFoam，我想着用sprayFoam应该也可以解决这个问题。 今天就花了很久在看sprayFoam 但是发现了一个问题，我可以把所有的化学反应关闭，但是在设置sprayCloudProperties的时候遇见了问题。我是基于 tutorials/lagrangian/sprayFoam/aachenBomb 这个case直接改的。这个case是喷入C7H16，我这里如果想只喷入固体颗粒，需要在singlePhaseMixtureCoeffs 下面的phases 里面写成 solid就可以了，然后再在 thermophysicalProperties 里面也改。但是运行的时候出现这个的时候出现 segmentation fault。

    请问什么问题。

  • chpjz0391

    我之前把例子中喷入的液体换成水倒是没什么问题？但是就是纯固体不行

  • 李东岳

    @chpjz0391 说：

    前把例子中喷入的液体换成水倒是没什么问题？但

    Hi，

    好图！非常直观

    或许把 segmentation fault这部分的log贴出来更方便诊断 :)

  • 散漫守望2016
    离散相副教授

    之前写的一个pimpleCoupledKinematicParcelFoam求解器，给你共享一下，大致下面这样，这个solver的主要部分，可以求解你的那个问题。
    主程序：
    int main(int argc, char *argv[])
    {
    argList::addOption
    (
    "cloudName",
    "name",
    "specify alternative cloud name. default is 'kinematicCloud'"
    );
    #include "setRootCase.H"
    #include "createTime.H"
    #include "createMesh.H"
    #include "readGravitationalAcceleration.H"
    #include "readTransportProperties.H"
    #include "createFields.H"
    #include "initContinuityErrs.H"
    #include "readTimeControls.H"
    #include "CourantNo.H"
    #include "setInitialDeltaT.H"

    pimpleControl pimple(mesh);

// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //

Info<< "\nStarting time loop\n" << endl;

while (runTime.run())
{

    #include "readTimeControls.H"
    #include "CourantNo.H"
    #include "setDeltaT.H"

    runTime++;

    Info<< "Time = " << runTime.timeName() << nl << endl;

    Info<< "Evolving " << kinematicCloud.name() << endl;
    kinematicCloud.evolve();
    kinematicCloud.info();

    // --- Pressure-velocity PIMPLE corrector loop
    while (pimple.loop())
    {
        #include "UEqn.H"

        // --- Pressure corrector loop
        while (pimple.correct())
        {
            #include "pEqn.H"
        }

        if (pimple.turbCorr())
        {
            turbulence->correct();
        }
    }
    
    runTime.write();

    Info<< "ExecutionTime = " << runTime.elapsedCpuTime() << " s"
        << "  ClockTime = " << runTime.elapsedClockTime() << " s"
        << nl << endl;
}

Info<< "End\n" << endl;

return 0;

    }

    CreateFields.H文件
    Info<< "Reading normalised dynamic presure field p_rgh\n" << endl;
    volScalarField p_rgh
    (
    IOobject
    (
    "p_rgh",
    runTime.timeName(),
    mesh,
    IOobject::MUST_READ,
    IOobject::AUTO_WRITE
    ),
    mesh
    );

    Info<< "Reading field U\n" << endl;
volVectorField U
(
    IOobject
    (
        "U",
        runTime.timeName(),
        mesh,
        IOobject::MUST_READ,
        IOobject::AUTO_WRITE
    ),
    mesh
);

Info<< "Creating phi [m3/s] " << endl;
#include "createPhi.H"

Info<< "Creating turbulence model\n" << endl;
singlePhaseTransportModel laminarTransport(U, phi);
autoPtr<incompressible::turbulenceModel> turbulence
(
    incompressible::turbulenceModel::New(U, phi, laminarTransport)
);

//Info<< "Creating source list\n" << endl;
fv::IOoptionList fvOptions(mesh);


Info<< "Calculating field g.h\n" << endl;
volScalarField gh("gh", g & mesh.C());
surfaceScalarField ghf("ghf", g & mesh.Cf());

Info<< "Creating absolute pressure field p [Pa]\n" << endl;
volScalarField p
(
    IOobject
    (
        "p",
        runTime.timeName(),
        mesh,
        IOobject::NO_READ,
        IOobject::AUTO_WRITE
    ),
    rhoVal*(p_rgh + gh)
);

Info<< "Reading field U\n" << endl;

label pRefCell = 0;
scalar pRefValue = 0.0;
setRefCell
(
    p,
    p_rgh,
    mesh.solutionDict().subDict("PIMPLE"),
    pRefCell,
    pRefValue
);

if (p_rgh.needReference())
{
    p += dimensionedScalar
    (
        "p",
        p.dimensions(),
        pRefValue - getRefCellValue(p, pRefCell)
    );
}

Info<< "Creating mass flux density surface field rhoPhi [kg/s] " << endl;
surfaceScalarField rhoPhi
(
    IOobject
    (
        "rhoPhi",
        runTime.timeName(),
        mesh,
        IOobject::READ_IF_PRESENT,
        IOobject::AUTO_WRITE
    ),
    rhoVal*linearInterpolate(U) & mesh.Sf()
);

// dummy-fields to satisfy the requirements of the kinematicCloud
Info<< "Creating rho field " << endl;
volScalarField rho
(
    IOobject
    (
        "rho",
        runTime.timeName(),
        mesh,
        IOobject::NO_READ,
        IOobject::NO_WRITE
    ),
    mesh,
    rhoVal,
    "zeroGradient"
);

Info<< "Creating mu field " << endl;
volScalarField mu
(
    IOobject
    (
        "mu",
        runTime.timeName(),
        mesh,
        IOobject::NO_READ,
        IOobject::NO_WRITE
    ),
    mesh,
    rhoVal*nu,
    "zeroGradient"
);

Info<< "Constructing kinematicCloud" << endl;
basicKinematicCollidingCloud kinematicCloud
(
    "kinematicCloud",
    rho,
    U,
    mu,
    g
);

    **
    readTransportProperties文件**：
    Info<< "\nReading transportProperties\n" << endl;

    IOdictionary transportProperties
(
    IOobject
    (
        "transportProperties",
        runTime.constant(),
        mesh,
        IOobject::MUST_READ_IF_MODIFIED,
        IOobject::NO_WRITE
    )
);

dimensionedScalar nu
(
    transportProperties.lookup("nu")
);

dimensionedScalar rhoVal
(
    transportProperties.lookup("rho")
);
Info<< "rhoVal is " << rhoVal << endl;

dimensionedScalar irhoVal("irhoVal",(1.0/rhoVal));

    **
    UEqn.H文件：**
    // Solve the momentum equation
    fvVectorMatrix UEqn
    (
    fvm::ddt(U)

    • fvm::div(phi, U)
    • turbulence->divDevReff(U)
      ==
      irhoVal*kinematicCloud.SU(U)
    • fvOptions(U)
      );

    UEqn.relax();

    fvOptions.constrain(UEqn);

    if (pimple.momentumPredictor())
    {
    solve
    (
    UEqn
    ==
    fvc::reconstruct
    (
    -fvc::snGrad(p_rgh)*mesh.magSf()
    )
    );
    }

    pEqn.H文件：
    {
    volScalarField rAU("rAU", 1.0/UEqn.A());
    surfaceScalarField rAUf("(1|A(U))", fvc::interpolate(rAU));

    U = rAU*UEqn.H();

    phi = (fvc::interpolate(U) & mesh.Sf());
    //+ fvc::ddtCorr(U, phi);

    adjustPhi(phi, U, p_rgh);

    while (pimple.correctNonOrthogonal())
    {
    fvScalarMatrix p_rghEqn
    (
    fvm::laplacian(rAU, p_rgh) == fvc::div(phi)
    );

      p_rghEqn.setReference(pRefCell, getRefCellValue(p_rgh, pRefCell));

  p_rghEqn.solve(mesh.solver(p_rgh.select(pimple.finalInnerIter())));

  if (pimple.finalNonOrthogonalIter())
  {
      // Calculate the conservative fluxes
      phi -= p_rghEqn.flux();

      // Explicitly relax pressure for momentum corrector
      p_rgh.relax();

      // Correct the momentum source with the pressure gradient flux
      // calculated from the relaxed pressure
      U -= rAU*fvc::reconstruct((p_rghEqn.flux())/rAUf);
      U.correctBoundaryConditions();
  }

    }

    #include "continuityErrs.H"

    p = rhoVal*(p_rgh + gh);

    if (p_rgh.needReference())
    {
    p += dimensionedScalar
    (
    "p",
    p.dimensions(),
    pRefValue - getRefCellValue(p, pRefCell)
    );
    p_rgh = irhoVal*p - gh;
    }
    }

  • chpjz0391

    这几天放假没在学校，才看到回复。谢谢分享 。我研究一下。万分感谢

Log in to reply