CFD Online Discussion Forums - Finding moving equation of parcel in icoLagrangianFoam

CFD Online Discussion Forums (https://www.cfd-online.com/Forums/)

- OpenFOAM Programming & Development (https://www.cfd-online.com/Forums/openfoam-programming-development/)

- - Finding moving equation of parcel in icoLagrangianFoam (https://www.cfd-online.com/Forums/openfoam-programming-development/241968-finding-moving-equation-parcel-icolagrangianfoam.html)

Finding moving equation of parcel in icoLagrangianFoam

Hi everyone,

I am a new foamer and I have a problem when adding electrical force in icoLagrangianFoam. It seems that there is no relation between void Foam::Cloud<ParticleType>::move(TrackingData& td) described in kinematicCloud.evolve() with moving equations of parcel.
Could anyone help me find the class that contains equations of parcel.

Thank you very much!

Hi Truong.

Particle equations of motion are found in $FOAM_SRC/lagrangian/intermediate/parcels/Templates/KinematicParcel/KinematicParcel.C (at least for OpenFOAM 6 anyway). Then in there, there is a function Foam::KinematicParcel<ParcelType>::calcVelocity. This function calculates the velocity increment based on the explicit and implicit forces acting on the particle:

Code:

    // Shortcut splitting assuming no implicit non-coupled force ...

    // Calculate the integration coefficients

    const vector acp = (Fcp.Sp()*td.Uc() + Fcp.Su())/massEff;

    const vector ancp = (Fncp.Su() + Su)/massEff;

    const scalar bcp = Fcp.Sp()/massEff;

Then the velocity increment is computed

Code:

    // Integrate to find the new parcel velocity

    const vector deltaU = cloud.UIntegrator().delta(U_, dt, acp + ancp, bcp);

    const vector deltaUncp = ancp*dt;

    const vector deltaUcp = deltaU - deltaUncp;



    // Calculate the new velocity and the momentum transfer terms

    vector Unew = U_ + deltaU;

There is a lot of templating, so it can take a while of following code to find where various properties are actually computed.

Quote:

Originally Posted by joshwilliams (Post 825177)

Particle equations of motion are found in $FOAM_SRC/lagrangian/intermediate/parcels/Templates/KinematicParcel/KinematicParcel.C (at least for OpenFOAM 6 anyway).

I should add that this is of course based on whichever force submodels you apply to the particle. For example, drag, lift, gravity etc found in

Code:

$FOAM_SRC/lagrangian/intermediate/submodels/Kinematic/ParticleForces

So for your electrical force, you should probably begin by copying one of the classes in this directory that best matches what you want to, and modify it from there.

Dear Josh Williams,

Thank you very much for your reply. You have provided me with valuable information.
I am a new foamer so it is helpful for me If you guide me how the function calcVelocity() in KinematicParcel.C related to function move(TrackingData& td) in CloudTemplate.C

Code:

void Foam::Cloud<ParticleType>::move(TrackingData& td)

{

    addProfile2(moveProfile,"Cloud<ParticleType>::move_"+this->name());



    const globalMeshData& pData = polyMesh_.globalData();

    const labelList& processorPatches = pData.processorPatches();

    const labelList& processorPatchIndices = pData.processorPatchIndices();

    const labelList& processorPatchNeighbours =

        pData.processorPatchNeighbours();



    // Initialise the setpFraction moved for the particles

    forAllIter(typename Cloud<ParticleType>, *this, pIter)

    {

        pIter().stepFraction() = 0;

    }



    // Assume there will be particles to transfer

    bool transfered = true;



    // While there are particles to transfer

    while (transfered)

    {

        // List of lists of particles to be transfered for all the processor

        // patches

        List<IDLList<ParticleType> > transferList(processorPatches.size());



        // Loop over all particles

        forAllIter(typename Cloud<ParticleType>, *this, pIter)

        {

            ParticleType& p = pIter();



            // Move the particle

            bool keepParticle = p.move(td);



            // If the particle is to be kept

            // (i.e. it hasn't passed through an inlet or outlet)

            if (keepParticle)

            {

                // If we are running in parallel and the particle is on a

                // boundary face

                if (Pstream::parRun() && p.facei_ >= pMesh().nInternalFaces())

                {

                    label patchi = pMesh().boundaryMesh().whichPatch(p.facei_);

                    label n = processorPatchIndices[patchi];



                    // ... and the face is on a processor patch

                    // prepare it for transfer

                    if (n != -1)

                    {

                        p.prepareForParallelTransfer(patchi, td);

                        transferList[n].append(this->remove(&p));

                    }

                }

            }

            else

            {

                deleteParticle(p);

            }

        }



        if (Pstream::parRun())

        {

            \\ something here

        }

        else

        {

            transfered = false;

        }

    }

}

This function is contained in function "kinematicCloud.evolve()" of icoLagrangianFoam.C

Code:

Info<< "Evolving kinematicCloud" << endl;

        kinematicCloud.evolve();

        kinematicCloud.info();



        fvVectorMatrix UEqn

        (

            fvm::ddt(U)

          + fvm::div(phi, U)

          - fvm::laplacian(nu, U)

            ==

            kinematicCloud.SU()/rho

        );



        solve(UEqn == -fvc::grad(p));

I think the move() function in "kinematicCloud.evolve()" is used to calculate the position of parcel, but I cannot find this relation between move() function and calcVelocity() function.

I hope that you could clarify this problem for me.
Thank you very much!

Hi Truong,

I remember trying to find this myself a year ago roughly. I think the relationship between kinematicCloud and kinematicParcel classes are very heavily templated and abstracted so it is a bit tricky to understand. But, I would say if you are simply adding an electrical force, you do not need to worry about this. You just simply copy a particle force template and modify it for your electrical force, no?

Best,
Josh

Hi Josh,

Thank you very much for your reply.
I got the link between KinematicCloud and KinematicPacel as you instructed me.

I have another problem and it is very kind if you could help me solve it.
I want to modify the calcVelocity function in KinematicParcel.C because this function controls the particle motion in move() function of CloudTemplate.C.
If you could remember this problem, please give me details for solving it.

Thank you very much once more time!

Best,
Truong