DPMFoam - Serious Error --particle-laden flow in simple geometric config

benz25 · February 18, 2015, 05:51

Dear Foamers,

I am currently studying multiphase flows in very small geometric configurations using OpenFOAM. I am quite new with OpenFOAM and thought of tryinng DPMFoam

.

My goal is to simulate particle-laden laminar flow in the following geometry (G5.png). I have designed the geometry and mesh in ANSYS Fluent, after which I converted the case file to OpenFOAM format and then ran numerical simulations in OpenFOAM.

The properties of the computational domain and particles are: 3 mm by 1 mm. Inlet velocity at 1 mm/s. Fluid: Air. Laminar Flow. Solid Particles are Polystyrene.

Moreover, I've modified the 0, constant, & system files as follows:

kinematicCloudProperties:

Code:

/*--------------------------------*- C++ -*----------------------------------*\
| =========                 |                                                 |
| \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox           |
|  \\    /   O peration     | Version:  2.3.x                                 |
|   \\  /    A nd           | Web:      www.OpenFOAM.org                      |
|    \\/     M anipulation  |                                                 |
\*---------------------------------------------------------------------------*/
FoamFile
{
    version     2.0;
    format      ascii;
    class       dictionary;
    location    "constant";
    object      particleProperties;
}
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //

solution
{
    active          true;
    coupled         true;
    transient       yes;
    cellValueSourceCorrection off;

    interpolationSchemes
    {
        rho.air         cell;
        U.air           cellPoint;
        mu.air          cell;
    }

    integrationSchemes
    {
        U               Euler;
    }

    sourceTerms
    {
        schemes
        {
            U semiImplicit 1;
        }
    }
}

constantProperties
{
    parcelTypeId 1;

    rhoMin          1e-15;
    minParticleMass 1e-15;

    rho0            1000;
    youngsModulus   3e9;
    poissonsRatio   0.34;

    constantVolume  false;

    alphaMax        0.99;
}

subModels
{
    particleForces
    {
        ErgunWenYuDrag
        {
            alphac alpha.air;
        }
        gravity;
    }

    injectionModels
    {
        model1
        {
            type            manualInjection;
            massTotal       0;
            parcelBasisType fixed;
            nParticle       500;
            SOI             0;
            positionsFile   "kinematicCloudPositions";
            U0              ( 0 0 0 );
            sizeDistribution
            {
                type        fixedValue;
                fixedValueDistribution
                {
                    value   0.00005;
                }
            }
        }
    }

    dispersionModel none;

    patchInteractionModel localInteraction;

    localInteractionCoeffs
    {
        patches
        (
            topwall
            {
                type rebound;
                e    0.97;
                mu   0.09;
            }
            bottomwall
            {
                type rebound;
                e    0.97;
                mu   0.09;
            }
            inlet
            {
                type rebound;
                e    0.97;
                mu   0.09;
            }
            outlet
            {
                type rebound;
                e    0.97;
                mu   0.09;
            }
            obstaclewalls
            {
                type rebound;
                e    0.97;
                mu   0.09;
            }
        );
    }

    StandardWallInteractionCoeffs
    {
        type rebound;
        e    0.97;
        mu   0.09;
    }

    heatTransferModel none;

    surfaceFilmModel none;

    collisionModel pairCollision;

    pairCollisionCoeffs
    {
        maxInteractionDistance  0.0001;

        writeReferredParticleCloud no;

        pairModel pairSpringSliderDashpot;

        pairSpringSliderDashpotCoeffs
        {
            useEquivalentSize   no;
            alpha               0.02;
            b                   1.5;
            mu                  0.10;
            cohesionEnergyDensity 0;
            collisionResolutionSteps 12;
        };

        wallModel wallSpringSliderDashpot;

        wallSpringSliderDashpotCoeffs
        {
            useEquivalentSize no;
            collisionResolutionSteps 12;
            youngsModulus   3e9;
            poissonsRatio   0.23;
            alpha           0.01;
            b               1.5;
            mu              0.09;
            cohesionEnergyDensity 0;
        };

        UName U.air;
    }

    stochasticCollisionModel none;

    radiation off;
}


cloudFunctions
{}


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

KinematicCloudPositions

Code:

/*--------------------------------*- C++ -*----------------------------------*\
| =========                 |                                                 |
| \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox           |
|  \\    /   O peration     | Version:  2.3.x                                 |
|   \\  /    A nd           | Web:      www.OpenFOAM.org                      |
|    \\/     M anipulation  |                                                 |
\*---------------------------------------------------------------------------*/
FoamFile
{
    version     2.0;
    format      ascii;
    class       vectorField;
    object      kinematicCloudPositions;
}
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //

(
(0.0005 0.0018 0)
)

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

g:

Code:

/*--------------------------------*- C++ -*----------------------------------*\
| =========                 |                                                 |
| \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox           |
|  \\    /   O peration     | Version:  2.3.x                                 |
|   \\  /    A nd           | Web:      www.OpenFOAM.org                      |
|    \\/     M anipulation  |                                                 |
\*---------------------------------------------------------------------------*/
FoamFile
{
    version     2.0;
    format      ascii;
    class       uniformDimensionedVectorField;
    location    "constant";
    object      g;
}
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //

dimensions      [0 1 -2 0 0 0 0];
value           (0 -9.81 0);


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

p:

Code:

/*--------------------------------*- C++ -*----------------------------------*\
| =========                 |                                                 |
| \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox           |
|  \\    /   O peration     | Version:  2.3.x                                 |
|   \\  /    A nd           | Web:      www.OpenFOAM.org                      |
|    \\/     M anipulation  |                                                 |
\*---------------------------------------------------------------------------*/
FoamFile
{
    version     2.0;
    format      ascii;
    class       volScalarField;
    object      p;
}
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //

dimensions      [0 2 -2 0 0 0 0];

internalField   uniform 0;

boundaryField
{
    obstaclewalls
    {
        type            fixedFluxPressure;
        phi             phi.air;
        value           $internalField;
    }

    bottomwall
    {
        type            fixedFluxPressure;
        phi             phi.air;
        value           $internalField;
    }

    topwall
    {
        type            fixedFluxPressure;
        phi             phi.air;
        value           $internalField;
    }

    outlet
    {
        type            fixedFluxPressure;
        phi             phi.air;
        value           $internalField;
    }

    inlet
    {
        type            fixedValue;
        phi             phi.air;
        value           $internalField;
    }

    frontAndBack
    {
        type            empty;
    }
}

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

U.air:

Code:

/*--------------------------------*- C++ -*----------------------------------*\
| =========                 |                                                 |
| \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox           |
|  \\    /   O peration     | Version:  2.3.x                                 |
|   \\  /    A nd           | Web:      www.OpenFOAM.org                      |
|    \\/     M anipulation  |                                                 |
\*---------------------------------------------------------------------------*/
FoamFile
{
    version     2.0;
    format      binary;
    class       volVectorField;
    location    "0";
    object      U.air;
}
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //

dimensions      [0 1 -1 0 0 0 0];

internalField   uniform (0 0 0);

boundaryField
{
    obstaclewalls
    {
        type            fixedValue;
        value           uniform (0 0 0);
    }

    bottomwall
    {
        type            fixedValue;
        value           uniform (0 0 0);
    }

    topwall
    {
        type            fixedValue;
        value           uniform (0 0 0);
    }

    inlet
    {
        type            pressureInletOutletVelocity;
        phi             phi.air;
        inletValue      uniform (0 0 0);
        value           uniform (0 0 0);
    }

    outlet
    {
        type            interstitialInletVelocity;
        inletVelocity   uniform (0.00015 0 0);
        value           uniform (0.00015 0 0);
        phi             phi.air;
        alpha           alpha.air;
    }

    frontAndBack
    {
        type            empty;
    }
}


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

controlDict:

Code:

/*--------------------------------*- C++ -*----------------------------------*\
| =========                 |                                                 |
| \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox           |
|  \\    /   O peration     | Version:  2.3.x                                 |
|   \\  /    A nd           | Web:      www.OpenFOAM.org                      |
|    \\/     M anipulation  |                                                 |
\*---------------------------------------------------------------------------*/
FoamFile
{
    version     2.0;
    format      ascii;
    class       dictionary;
    location    "system";
    object      controlDict;
}
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //

application     DPMFoam;

startFrom       latestTime;

startTime       0;

stopAt          endTime;

endTime         5;

deltaT          2e-5;

writeControl    runTime;

writeInterval   0.01;

purgeWrite      0;

writeFormat     binary;

writePrecision  6;

writeCompression uncompressed;

timeFormat      general;

timePrecision   6;

runTimeModifiable yes;

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

fvSchemes:

Code:

/*--------------------------------*- C++ -*----------------------------------*\
| =========                 |                                                 |
| \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox           |
|  \\    /   O peration     | Version:  2.3.x                                 |
|   \\  /    A nd           | Web:      www.OpenFOAM.org                      |
|    \\/     M anipulation  |                                                 |
\*---------------------------------------------------------------------------*/
FoamFile
{
    version     2.0;
    format      ascii;
    class       dictionary;
    object      fvSchemes;
}
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //

ddtSchemes
{
    default Euler;
}

gradSchemes
{
    default         Gauss linear;
}

divSchemes
{
    default         none;

    div(alphaPhic,U.air)    Gauss linearUpwindV unlimited;
    div(((alpha.air*nuEff.air)*dev2(T(grad(U.air))))) Gauss linear;
}

laplacianSchemes
{
    default         Gauss linear corrected;
}

interpolationSchemes
{
    default         linear;
}

snGradSchemes
{
    default         corrected;
}

fluxRequired
{
    default         no;
    p;
    kinematicCloud:theta;
}

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

fvSolution

Code:

/*--------------------------------*- C++ -*----------------------------------*\
| =========                 |                                                 |
| \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox           |
|  \\    /   O peration     | Version:  2.3.x                                 |
|   \\  /    A nd           | Web:      www.OpenFOAM.org                      |
|    \\/     M anipulation  |                                                 |
\*---------------------------------------------------------------------------*/
FoamFile
{
    version     2.0;
    format      ascii;
    class       dictionary;
    location    "system";
    object      fvSolution;
}
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //

solvers
{
    "(p|kinematicCloud:theta)"
    {
        solver          GAMG;
        tolerance       1e-06;
        relTol          0.01;
        smoother        GaussSeidel;
        cacheAgglomeration true;
        nCellsInCoarsestLevel 10;
        agglomerator    faceAreaPair;
        mergeLevels     1;
    }

    "(p|kinematicCloud:theta)Final"
    {
        solver          GAMG;
        tolerance       1e-06;
        relTol          0;
        smoother        GaussSeidel;
        cacheAgglomeration true;
        nCellsInCoarsestLevel 10;
        agglomerator    faceAreaPair;
        mergeLevels     1;
    }

    "(U.air|k|omega)"
    {
        solver          smoothSolver;
        smoother        symGaussSeidel;
        tolerance       1e-05;
        relTol          0.1;
    }

    "(U.air|k|omega)Final"
    {
        solver          smoothSolver;
        smoother        symGaussSeidel;
        tolerance       1e-05;
        relTol          0;
    }
}

PIMPLE
{
    nOuterCorrectors 1;
    nCorrectors     2;
    momentumPredictor yes;
    nNonOrthogonalCorrectors 0;
    pRefCell        0;
    pRefValue       0;
}

relaxationFactors
{
}


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

After running DPMFoam I get the following error...the solution diverges and the time-step continuity errors/Courant number gets extremely large after every iteration executed, as well as the linear momentum and energies! I highlighted these errors in red.

I initially thought I input the incorrect particle properties but that wasn't the case. So I suspect it's to do with the properties in KinematicCloud file. I 'm not sure if I've input the correct details in "injectionModels". For instance, since I've designed the geometry in OpenFOAM with the dimensions above, I'm not sure if it's to do with the incorrect coordinates for point of injection. I'm not sure where the error lies and I would be extremely grateful if someone could guide me to right direction..

Many thanks in advance

Regards,

benz25 · February 18, 2015, 05:57

This is the Error I got when running DPMFoam:

Code:

 

 .............

 .............

    From function solution::solver(const word&)
    in file matrices/solution/solution.C at line 381
    Lookup solver for pFinal
GAMG:  Solving for p, Initial residual = 0.00803955, Final residual = 5.66417e-07, No Iterations 18
Lookup fluxRequired for p
Field relaxation factor for pFinal is unset
Lookup fluxRequired for p
Lookup ddtScheme for ddt(alpha.air)
time step continuity errors : sum local = 0.201158, global = -0.0096157, cumulative = 0.258732
ExecutionTime = 79.97 s  ClockTime = 86 s

Courant Number mean: 11341.4 max: 658250
Time = 0.00056

Evolving kinematicCloud
Cache: Deleting volPointInterpolate(U.air), originating from U.air event No. 22417
Cache: Recalculating volPointInterpolate(U.air), originating from U.air event No. 22417
Cache: Storing volPointInterpolate(U.air), originating from U.air event No. 22417

Solving 2-D cloud kinematicCloud
    483 move-collide subCycles
Cloud: kinematicCloud
    Current number of parcels       = 1
    Current mass in system          = 1.65326e-10
    Linear momentum                 = (4.20393e-08 -5.42428e-09 0)
   |Linear momentum|                = 4.23878e-08
    Linear kinetic energy           = 5.43388e-06
    model1:
        number of parcels added     = 1
        mass introduced             = 1.65326e-10
    Parcels absorbed into film      = 0
    New film detached parcels       = 0
    Parcel fate (number, mass)      : patch topwall
      - escape                      = 0, 0
      - stick                       = 0, 0
    Parcel fate (number, mass)      : patch bottomwall
      - escape                      = 0, 0
      - stick                       = 0, 0
    Parcel fate (number, mass)      : patch inlet
      - escape                      = 0, 0
      - stick                       = 0, 0
    Parcel fate (number, mass)      : patch outlet
      - escape                      = 0, 0
      - stick                       = 0, 0
    Parcel fate (number, mass)      : patch obstaclewalls
      - escape                      = 0, 0
      - stick                       = 0, 0
    Rotational kinetic energy       = 9.52595e-10

Lookup interpolationScheme for interpolate(alpha.air)
Lookup gradScheme for grad(U.air)
Cache: Calculating grad(U.air), originating from U.air event No. 22417
Lookup divScheme for div(((alpha.air*nuEff.air)*dev2(T(grad(U.air)))))
Lookup laplacianScheme for laplacian((alpha.air*nuEff.air),U.air)
Lookup fluxRequired for U.air
Lookup gradScheme for grad(U.air)
Cache: Calculating grad(U.air), originating from U.air event No. 22417
Lookup ddtScheme for ddt(alpha.air)
Lookup divScheme for div(alphaPhic,U.air)
Lookup gradScheme for unlimited
Lookup ddtScheme for default
Cache: Calculating unlimited, originating from U.air event No. 22417
Lookup ddtScheme for ddt(alpha.air,U.air)
Find equation relaxation factor for U.airFinal
Lookup interpolationScheme for interpolate((1|A(U.air)))
Lookup interpolationScheme for interpolate((((1|A(U.air))*cloudVolSUSu)|rho.air))
Lookup snGradScheme for snGrad(p)
Lookup gradScheme for grad(p)
Cache: Calculating grad(p), originating from p event No. 22344

    From function solution::solverDict(const word&)
    in file matrices/solution/solution.C at line 369
    Lookup solver for U.airFinal
smoothSolver:  Solving for U.airx, Initial residual = 0.607062, Final residual = 0.00924022, No Iterations 1000
smoothSolver:  Solving for U.airy, Initial residual = 0.700008, Final residual = 0.00198385, No Iterations 1000
Lookup ddtScheme for ddt(U.air)
Lookup interpolationScheme for interpolate(U.air_0)
Lookup interpolationScheme for interpolate(HbyA)
Lookup ddtScheme for ddt(alpha.air)
Lookup laplacianScheme for laplacian((alphacf*Dp),p)
Lookup fluxRequired for p
Lookup gradScheme for grad(p)
Cache: Calculating grad(p), originating from p event No. 22864

    From function solution::solver(const word&)
    in file matrices/solution/solution.C at line 381
    Lookup solver for p
GAMG:  Solving for p, Initial residual = 0.889519, Final residual = 0.00496431, No Iterations 11
Lookup fluxRequired for p
Field relaxation factor for pFinal is unset
Lookup fluxRequired for p
Lookup ddtScheme for ddt(alpha.air)
time step continuity errors : sum local = 122.507, global = -5.7043, cumulative = -5.44557
Lookup ddtScheme for ddt(U.air)
Lookup interpolationScheme for interpolate(U.air_0)
Lookup interpolationScheme for interpolate(HbyA)
Lookup ddtScheme for ddt(alpha.air)
Lookup laplacianScheme for laplacian((alphacf*Dp),p)
Lookup fluxRequired for p
Lookup gradScheme for grad(p)
Cache: Calculating grad(p), originating from p event No. 23102

    From function solution::solver(const word&)
    in file matrices/solution/solution.C at line 381
    Lookup solver for pFinal
GAMG:  Solving for p, Initial residual = 0.0204461, Final residual = 6.77139e-07, No Iterations 19
Lookup fluxRequired for p
Field relaxation factor for pFinal is unset
Lookup fluxRequired for p
Lookup ddtScheme for ddt(alpha.air)
time step continuity errors : sum local = 3.27523, global = 0.143837, cumulative = -5.30173
ExecutionTime = 84.84 s  ClockTime = 91 s

Courant Number mean: 400882 max: 1.36343e+07
Time = 0.00058

Evolving kinematicCloud
Cache: Deleting volPointInterpolate(U.air), originating from U.air event No. 23242
Cache: Recalculating volPointInterpolate(U.air), originating from U.air event No. 23242
Cache: Storing volPointInterpolate(U.air), originating from U.air event No. 23242

Solving 2-D cloud kinematicCloud
    682 move-collide subCycles
Cloud: kinematicCloud
    Current number of parcels       = 1
    Current mass in system          = 1.65326e-10
    Linear momentum                 = (2.08372e-06 -1.09387e-06 0)
   |Linear momentum|                = 2.35339e-06
    Linear kinetic energy           = 0.01675
    model1:
        number of parcels added     = 1
        mass introduced             = 1.65326e-10
    Parcels absorbed into film      = 0
    New film detached parcels       = 0
    Parcel fate (number, mass)      : patch topwall
      - escape                      = 0, 0
      - stick                       = 0, 0
    Parcel fate (number, mass)      : patch bottomwall
      - escape                      = 0, 0
      - stick                       = 0, 0
    Parcel fate (number, mass)      : patch inlet
      - escape                      = 0, 0
      - stick                       = 0, 0
    Parcel fate (number, mass)      : patch outlet
      - escape                      = 0, 0
      - stick                       = 0, 0
    Parcel fate (number, mass)      : patch obstaclewalls
      - escape                      = 0, 0
      - stick                       = 0, 0
    Rotational kinetic energy       = 4.18539e-10

Lookup interpolationScheme for interpolate(alpha.air)
Lookup gradScheme for grad(U.air)
Cache: Calculating grad(U.air), originating from U.air event No. 23242
Lookup divScheme for div(((alpha.air*nuEff.air)*dev2(T(grad(U.air)))))
Lookup laplacianScheme for laplacian((alpha.air*nuEff.air),U.air)
Lookup fluxRequired for U.air
Lookup gradScheme for grad(U.air)
Cache: Calculating grad(U.air), originating from U.air event No. 23242
Lookup ddtScheme for ddt(alpha.air)
Lookup divScheme for div(alphaPhic,U.air)
Lookup gradScheme for unlimited
Lookup ddtScheme for default
Cache: Calculating unlimited, originating from U.air event No. 23242
Lookup ddtScheme for ddt(alpha.air,U.air)
Find equation relaxation factor for U.airFinal
Lookup interpolationScheme for interpolate((1|A(U.air)))
Lookup interpolationScheme for interpolate((((1|A(U.air))*cloudVolSUSu)|rho.air))
Lookup snGradScheme for snGrad(p)
Lookup gradScheme for grad(p)
Cache: Calculating grad(p), originating from p event No. 23169

    From function solution::solverDict(const word&)
    in file matrices/solution/solution.C at line 369
    Lookup solver for U.airFinal
smoothSolver:  Solving for U.airx, Initial residual = 0.56959, Final residual = 2.43449, No Iterations 1000
smoothSolver:  Solving for U.airy, Initial residual = 0.741332, Final residual = 9.96388e-06, No Iterations 994
Lookup ddtScheme for ddt(U.air)
Lookup interpolationScheme for interpolate(U.air_0)
Lookup interpolationScheme for interpolate(HbyA)
Lookup ddtScheme for ddt(alpha.air)
Lookup laplacianScheme for laplacian((alphacf*Dp),p)
Lookup fluxRequired for p
Lookup gradScheme for grad(p)
Cache: Calculating grad(p), originating from p event No. 23689

    From function solution::solver(const word&)
    in file matrices/solution/solution.C at line 381
    Lookup solver for p
GAMG:  Solving for p, Initial residual = 0.999993, Final residual = 0.008442, No Iterations 12
Lookup fluxRequired for p
Field relaxation factor for pFinal is unset
Lookup fluxRequired for p
Lookup ddtScheme for ddt(alpha.air)
time step continuity errors : sum local = 7.8273e+06, global = -254411, cumulative = -254417
Lookup ddtScheme for ddt(U.air)
Lookup interpolationScheme for interpolate(U.air_0)
Lookup interpolationScheme for interpolate(HbyA)
Lookup ddtScheme for ddt(alpha.air)
Lookup laplacianScheme for laplacian((alphacf*Dp),p)
Lookup fluxRequired for p
Lookup gradScheme for grad(p)
Cache: Calculating grad(p), originating from p event No. 23927

    From function solution::solver(const word&)
    in file matrices/solution/solution.C at line 381
    Lookup solver for pFinal
GAMG:  Solving for p, Initial residual = 0.0702429, Final residual = 6.4113e-07, No Iterations 22
Lookup fluxRequired for p
Field relaxation factor for pFinal is unset
Lookup fluxRequired for p
Lookup ddtScheme for ddt(alpha.air)
time step continuity errors : sum local = 45069.2, global = 1363.5, cumulative = -253053
ExecutionTime = 89.77 s  ClockTime = 96 s

Courant Number mean: 2.13843e+10 max: 5.23752e+11
Time = 0.0006

Evolving kinematicCloud
Cache: Deleting volPointInterpolate(U.air), originating from U.air event No. 24067
Cache: Recalculating volPointInterpolate(U.air), originating from U.air event No. 24067
Cache: Storing volPointInterpolate(U.air), originating from U.air event No. 24067

Solving 2-D cloud kinematicCloud
    1517 move-collide subCycles
Cloud: kinematicCloud
    Current number of parcels       = 1
    Current mass in system          = 1.65326e-10
    Linear momentum                 = (-7.94201 -1.06758 0)
   |Linear momentum|                = 8.01344
    Linear kinetic energy           = 1.94207e+11
    model1:
        number of parcels added     = 1
        mass introduced             = 1.65326e-10
    Parcels absorbed into film      = 0
    New film detached parcels       = 0
    Parcel fate (number, mass)      : patch topwall
      - escape                      = 0, 0
      - stick                       = 0, 0
    Parcel fate (number, mass)      : patch bottomwall
      - escape                      = 0, 0
      - stick                       = 0, 0
    Parcel fate (number, mass)      : patch inlet
      - escape                      = 0, 0
      - stick                       = 0, 0
    Parcel fate (number, mass)      : patch outlet
      - escape                      = 0, 0
      - stick                       = 0, 0
    Parcel fate (number, mass)      : patch obstaclewalls
      - escape                      = 0, 0
      - stick                       = 0, 0
    Rotational kinetic energy       = 785.921

Lookup interpolationScheme for interpolate(alpha.air)
Lookup gradScheme for grad(U.air)
Cache: Calculating grad(U.air), originating from U.air event No. 24067
Lookup divScheme for div(((alpha.air*nuEff.air)*dev2(T(grad(U.air)))))
Lookup laplacianScheme for laplacian((alpha.air*nuEff.air),U.air)
Lookup fluxRequired for U.air
Lookup gradScheme for grad(U.air)
Cache: Calculating grad(U.air), originating from U.air event No. 24067
Lookup ddtScheme for ddt(alpha.air)
Lookup divScheme for div(alphaPhic,U.air)
Lookup gradScheme for unlimited
Lookup ddtScheme for default
Cache: Calculating unlimited, originating from U.air event No. 24067
Lookup ddtScheme for ddt(alpha.air,U.air)
Find equation relaxation factor for U.airFinal
Lookup interpolationScheme for interpolate((1|A(U.air)))
Lookup interpolationScheme for interpolate((((1|A(U.air))*cloudVolSUSu)|rho.air))
Lookup snGradScheme for snGrad(p)
Lookup gradScheme for grad(p)
Cache: Calculating grad(p), originating from p event No. 23994

    From function solution::solverDict(const word&)
    in file matrices/solution/solution.C at line 369
    Lookup solver for U.airFinal
smoothSolver:  Solving for U.airx, Initial residual = 0.609695, Final residual = 15522.5, No Iterations 1000
smoothSolver:  Solving for U.airy, Initial residual = 0.654515, Final residual = 2.97979e-05, No Iterations 1000
Lookup ddtScheme for ddt(U.air)
Lookup interpolationScheme for interpolate(U.air_0)
Lookup interpolationScheme for interpolate(HbyA)
Lookup ddtScheme for ddt(alpha.air)
Lookup laplacianScheme for laplacian((alphacf*Dp),p)
Lookup fluxRequired for p
Lookup gradScheme for grad(p)
Cache: Calculating grad(p), originating from p event No. 24514

    From function solution::solver(const word&)
    in file matrices/solution/solution.C at line 381
    Lookup solver for p
GAMG:  Solving for p, Initial residual = 1, Final residual = 0.00838736, No Iterations 11
Lookup fluxRequired for p
Field relaxation factor for pFinal is unset
Lookup fluxRequired for p
Lookup ddtScheme for ddt(alpha.air)
time step continuity errors : sum local = 1.10465e+15, global = -2.77922e+13, cumulative = -2.77922e+13
Lookup ddtScheme for ddt(U.air)
Lookup interpolationScheme for interpolate(U.air_0)
Lookup interpolationScheme for interpolate(HbyA)
Lookup ddtScheme for ddt(alpha.air)
Lookup laplacianScheme for laplacian((alphacf*Dp),p)
Lookup fluxRequired for p
Lookup gradScheme for grad(p)
Cache: Calculating grad(p), originating from p event No. 24752

    From function solution::solver(const word&)
    in file matrices/solution/solution.C at line 381
    Lookup solver for pFinal
GAMG:  Solving for p, Initial residual = 0.0746153, Final residual = 4.38279e-07, No Iterations 24
Lookup fluxRequired for p
Field relaxation factor for pFinal is unset
Lookup fluxRequired for p
Lookup ddtScheme for ddt(alpha.air)
time step continuity errors : sum local = 3.38976e+12, global = 1.15446e+11, cumulative = -2.76767e+13
ExecutionTime = 135.44 s  ClockTime = 141 s

Courant Number mean: 2.27238e+18 max: 1.25889e+20
Time = 0.00062

Evolving kinematicCloud
Cache: Deleting volPointInterpolate(U.air), originating from U.air event No. 24892
Cache: Recalculating volPointInterpolate(U.air), originating from U.air event No. 24892
Cache: Storing volPointInterpolate(U.air), originating from U.air event No. 24892

Solving 2-D cloud kinematicCloud
    30702 move-collide subCycles
--> FOAM Warning :
    From function Time::operator++()
    in file db/Time/Time.C at line 1055
    Increased the timePrecision from 6 to 7 to distinguish between timeNames at time 0.000600001
Cloud: kinematicCloud
    Current number of parcels       = 1
    Current mass in system          = 1.65326e-10
    Linear momentum                 = (-5.25941e+08 -1.3121e+08 0)
   |Linear momentum|                = 5.4206e+08
    Linear kinetic energy           = 8.88635e+26
    model1:
        number of parcels added     = 1
        mass introduced             = 1.65326e-10
    Parcels absorbed into film      = 0
    New film detached parcels       = 0
    Parcel fate (number, mass)      : patch topwall
      - escape                      = 0, 0
      - stick                       = 0, 0
    Parcel fate (number, mass)      : patch bottomwall
      - escape                      = 0, 0
      - stick                       = 0, 0
    Parcel fate (number, mass)      : patch inlet
      - escape                      = 0, 0
      - stick                       = 0, 0
    Parcel fate (number, mass)      : patch outlet
      - escape                      = 0, 0
      - stick                       = 0, 0
    Parcel fate (number, mass)      : patch obstaclewalls
      - escape                      = 0, 0
      - stick                       = 0, 0
    Rotational kinetic energy       = 785.921

Lookup interpolationScheme for interpolate(alpha.air)
Lookup gradScheme for grad(U.air)
Cache: Calculating grad(U.air), originating from U.air event No. 24892
Lookup divScheme for div(((alpha.air*nuEff.air)*dev2(T(grad(U.air)))))
Lookup laplacianScheme for laplacian((alpha.air*nuEff.air),U.air)
Lookup fluxRequired for U.air
Lookup gradScheme for grad(U.air)
Cache: Calculating grad(U.air), originating from U.air event No. 24892
Lookup ddtScheme for ddt(alpha.air)
Lookup divScheme for div(alphaPhic,U.air)
Lookup gradScheme for unlimited
Lookup ddtScheme for default
Cache: Calculating unlimited, originating from U.air event No. 24892
Lookup ddtScheme for ddt(alpha.air,U.air)
Find equation relaxation factor for U.airFinal
Lookup interpolationScheme for interpolate((1|A(U.air)))
Lookup interpolationScheme for interpolate((((1|A(U.air))*cloudVolSUSu)|rho.air))
Lookup snGradScheme for snGrad(p)
Lookup gradScheme for grad(p)
Cache: Calculating grad(p), originating from p event No. 24819

    From function solution::solverDict(const word&)
    in file matrices/solution/solution.C at line 369
    Lookup solver for U.airFinal
#0  Foam::error::printStack(Foam::Ostream&) at ??:?
#1  Foam::sigFpe::sigHandler(int) at ??:?
#2   in "/lib64/libc.so.6"
#3  Foam::symGaussSeidelSmoother::smooth(Foam::word const&, Foam::Field<double>&, Foam::lduMatrix const&, Foam::Field<double> const&, Foam::FieldField<Foam::Field, double> const&, Foam::UPtrList<Foam::lduInterfaceField const> const&, unsigned char, int) at ??:?
#4  Foam::symGaussSeidelSmoother::smooth(Foam::Field<double>&, Foam::Field<double> const&, unsigned char, int) const at ??:?
#5  Foam::smoothSolver::solve(Foam::Field<double>&, Foam::Field<double> const&, unsigned char) const at ??:?
#6  Foam::fvMatrix<Foam::Vector<double> >::solveSegregated(Foam::dictionary const&) at ??:?
#7  Foam::fvMatrix<Foam::Vector<double> >::solve(Foam::dictionary const&) at ??:?
#8  Foam::SolverPerformance<double> Foam::solve<Foam::Vector<double> >(Foam::tmp<Foam::fvMatrix<Foam::Vector<double> > > const&) at ??:?
#9
 at ??:?
#10  __libc_start_main in "/lib64/libc.so.6"
#11
 at /usr/src/packages/BUILD/glibc-2.11.3/csu/../sysdeps/x86_64/elf/start.S:116
Floating point exception

ansubru · March 16, 2015, 05:10

Hi!!

From what i can see, you dont have any problems with particle ínjection (if you did, your simulation will crash at the beginning)... The reason for your crash is that some of the bc's are incorrectly defined (especially the p bc's) as indicated by the error flags..

Code:

From function solution::solverDict(const word&)
    in file matrices/solution/solution.C at line 369
    Lookup solver for U.airFinal
#0  Foam::error::printStack(Foam::Ostream&) at ??:?
#1  Foam::sigFpe::sigHandler(int) at ??:?
#2   in "/lib64/libc.so.6"
#3  Foam::symGaussSeidelSmoother::smooth(Foam::word const&, Foam::Field<double>&, Foam::lduMatrix const&, Foam::Field<double> const&, Foam::FieldField<Foam::Field, double> const&, Foam::UPtrList<Foam::lduInterfaceField const> const&, unsigned char, int) at ??:?
#4  Foam::symGaussSeidelSmoother::smooth(Foam::Field<double>&, Foam::Field<double> const&, unsigned char, int) const at ??:?
#5  Foam::smoothSolver::solve(Foam::Field<double>&, Foam::Field<double> const&, unsigned char) const at ??:?
#6  Foam::fvMatrix<Foam::Vector<double> >::solveSegregated(Foam::dictionary const&) at ??:?
#7  Foam::fvMatrix<Foam::Vector<double> >::solve(Foam::dictionary const&) at ??:?
#8  Foam::SolverPerformance<double> Foam::solve<Foam::Vector<double> >(Foam::tmp<Foam::fvMatrix<Foam::Vector<double> > > const&) at ??:?
#9
 at ??:?

You cannot define a 'pressureinletoutlet' or ínterstitialinletoutlet' boundary as a fixedfluxpressure boundary... in such case rename the corresponding p boundary as 'fixedValue'... this might solve your problem...

Additionally, i can see that your particle size is 0.00005m, hence correspondingly the maxinteractiondistance defined in the kinematiccloudproperties file must also be 0.00005.

Additionally, if you want better control over your simulations, i would advice you to setup your case from scratch,i.e geometry and meshing, in OF (your setup is relatively simple)... In this way you are sure that the solutions are not affected by the different ways in which your case has been meshed...

Hope this msg will help you

Kind regards

ansubru

benz25 · March 19, 2015, 09:14

Hi,

Thank you for your reply. Indeed the B.C.s have been incorrectly assigned. I changed the P inlet & outlet boundaries to fixedValue.

I also designed a simple clear channel geometry and mesh in OF (see attached file ee22.png). The same B.C.s and properties have been applied (got rid of bc's for obstacle wall) with corrections to he P inlet/outlet B.C.s as you've mentioned.

After executing DPMFoam, I get the error (see ERROR Injection.txt). Someone raised this query before: http://www.cfd-online.com/Forums/ope...cro-scale.html

He solved it by assigning an injection point away from the boundary planes.

I tried that but I still get the same error

Not too sure where I went wrong?? Any help greatly appreciated

Many thanks & BR,

Kevin

Quote:

Originally Posted by ansubru

Hi!!

From what i can see, you dont have any problems with particle ínjection (if you did, your simulation will crash at the beginning)... The reason for your crash is that some of the bc's are incorrectly defined (especially the p bc's) as indicated by the error flags..

Code:

From function solution::solverDict(const word&)
    in file matrices/solution/solution.C at line 369
    Lookup solver for U.airFinal
#0  Foam::error::printStack(Foam::Ostream&) at ??:?
#1  Foam::sigFpe::sigHandler(int) at ??:?
#2   in "/lib64/libc.so.6"
#3  Foam::symGaussSeidelSmoother::smooth(Foam::word const&, Foam::Field<double>&, Foam::lduMatrix const&, Foam::Field<double> const&, Foam::FieldField<Foam::Field, double> const&, Foam::UPtrList<Foam::lduInterfaceField const> const&, unsigned char, int) at ??:?
#4  Foam::symGaussSeidelSmoother::smooth(Foam::Field<double>&, Foam::Field<double> const&, unsigned char, int) const at ??:?
#5  Foam::smoothSolver::solve(Foam::Field<double>&, Foam::Field<double> const&, unsigned char) const at ??:?
#6  Foam::fvMatrix<Foam::Vector<double> >::solveSegregated(Foam::dictionary const&) at ??:?
#7  Foam::fvMatrix<Foam::Vector<double> >::solve(Foam::dictionary const&) at ??:?
#8  Foam::SolverPerformance<double> Foam::solve<Foam::Vector<double> >(Foam::tmp<Foam::fvMatrix<Foam::Vector<double> > > const&) at ??:?
#9
 at ??:?

You cannot define a 'pressureinletoutlet' or ínterstitialinletoutlet' boundary as a fixedfluxpressure boundary... in such case rename the corresponding p boundary as 'fixedValue'... this might solve your problem...

Additionally, i can see that your particle size is 0.00005m, hence correspondingly the maxinteractiondistance defined in the kinematiccloudproperties file must also be 0.00005.

Additionally, if you want better control over your simulations, i would advice you to setup your case from scratch,i.e geometry and meshing, in OF (your setup is relatively simple)... In this way you are sure that the solutions are not affected by the different ways in which your case has been meshed...

Hope this msg will help you

Kind regards

ansubru

ansubru · March 19, 2015, 09:37

Hi Kevin!!

I have a few suggestions for you :

--> Firstly run the command çheckMesh' to see mesh statistics... I have a feeling you have not followed the 'right hand rule' while assigning patches and that your mesh is quite poor (please read the blockMesh help in detail before you begin meshing)... refer here..http://www.openfoam.org/docs/user/blockMesh.php

--> Secondly, you need to define your injection correctly, if your mesh has been 'drawn' in a different co-ordinate system, you could have issues when injections have to be initialized... And try to use alternate particle injection strategies.. I use patchinjection, and that seems to work quite well for me

...

--> Thirdly, try not to use the exact young's modulus while computing... Use a value 1000 times lower (for example try not to exceed the value used in the DPMFoam tutorial)... Your simulation will still be correct, dont worry.. Its an approximation almost all DEM works use so as to ease the computational burden... I have commented on a post regarding the same.. read more here..

http://www.cfd-online.com/Forums/ope...tml#post513475

Hope that helped!!

Regards

ansubru

benz25 · March 25, 2015, 08:40

Dear Ananda,

Thank you very much for your assistance! I went through the websites you posted and took note of the meshing procedure in OF. By following the right hand rule, and using patchinjection, I was able to execute the solver without any issues

. I'm going to try cone injection model and try and modify particle properties using different size distribution model and see how it goes...

However, in the KinematicCloudProperties file, I'm still not sure what the following terms mean and I'm not sure how or under what situations should I modify those terms (highlighted in red). DPMFoam uses hard sphere model and Hertzian contact if i'm not mistaken?..According to http://www.cfd-online.com/Forums/ope...e-dpmfoam.html does this mean a non adhesive Hertzian contact model is used or adhesive contact i.e. JKR/DMT? Can we select other contact models?

And what is the difference between using the openfoam solver DPMFoam and CFDEM (liggghts+openfoam) when it comes to simulating particle-fluid flows? According to website, both codes have 4-way coupling functionality i.e. particle-particle and particle-fluid interactions. So i'm not sure if one of them has something the other doesn't have? Any distinct difference between them

?

Many thanks again for your assistance!!!

Kevin

Code:

pairCollisionCoeffs
    {
        maxInteractionDistance  0.000005;

        writeReferredParticleCloud no;

        pairModel pairSpringSliderDashpot;

        pairSpringSliderDashpotCoeffs
        {
            useEquivalentSize   no;
            alpha               0.02;
            b                   1.5;
            mu                  0.10;
            cohesionEnergyDensity 0;
            collisionResolutionSteps 12;
        };

        wallModel wallSpringSliderDashpot;

        wallSpringSliderDashpotCoeffs
        {
            useEquivalentSize no;
            collisionResolutionSteps 12;
            youngsModulus   1e4;
            poissonsRatio   0.23;
            alpha           0.01;
            b               1.5;
            mu              0.09;
            cohesionEnergyDensity 0;
        };

        UName U.air;
    }

ansubru · March 25, 2015, 11:19

Hi Kevin!!

Firstly, DPM-Foam is a 'soft sphere' contact model based solver. The entire concept of CFD-DEM has been developed by Tsuji et.al (look at the paper below) while studying plug flow conveying. All relevant parameters for the method of DEM are detailed in their paper, you can have a detailed read and I am sure that will answer all your questions...

https://www.dropbox.com/sh/6bis355hp...ihf45kyra?dl=0

Equivalent size is a concept that comes into play when you have non spherical particles. My suggestion would be for you to read about hertzian contact models and why they are universally accepted as one of the best contact models available for simulating multi-particle contact. Additionally, the DPMFoam solver can also account for cohension, but I am not very aware of the models behind that. The config file you have shared is for 'non-adhesive' contacts. As of now, DPMFoam only uses linear and hertzian contact pair spring slider dashpot models for contact. Any new model, will unfortunately have to be developed you...

And finally, about the differences between the various solvers (i.e CFDEM-LIGGGHTS, etc.) there is'nt much. They have the same features and are equally effective in resolving granular flow physics. LIGGGHTS has one major advantage that there are lot of supporting documentation/tutorials on setting up and analyzing a simulation. Apart from that, there is also the added capability of modelling dynamic meshes- moving geometries (which unfortunately cant be done with DPMFoam right now)... However, DPMFoam is more malleable to your needs than LIGGGHTS, i.e if you of course read and understand the source code. Else both tools are equally capable and its a matter of choice really. Personally, I couldnt get LIGGGHTS to work

some incompatibility with my hardware, hence I adopted DPMFoam

Hope that helped

All the best in your quest!!

ansubru

benz25 · April 2, 2015, 12:01

Hi Ananda,

Thank you very much for your reply. I did have some compatibility issues with liggghts as well. After spending some time with the installation, I did manage to install it using Ubuntu 12 LTS and OpenFOAM 2.1.1. I used either Helyx or Salome for meshing.

Thanks for sharing the document pertaining to CFD-DEM developed by Tsuji (1993)

. I've read two of Tsuji's papers. I see they use the spring dashpot model for the forces whilst the Hertz theory is used to evaluate the stiffness. So DPMFoam allows users to select either linear-spring–dashpot model or
non-linear Hertzian-spring–dashpot by setting b=1 or 1.5 respectively:

Code:

[OpenFOAM-2.3.x/src/lagrangian/intermediate/submodels/Kinematic/CollisionModel/PairCollision/PairModel/PairSpringSliderDashpot/PairSpringSliderDashpot.H

Interestingly, after going through additional papers regarding surface energies (used in cohesion models), and then looking back at the DPMFoam code, I noted that in Line 217

Code:

OpenFOAM-2.3.x/src/lagrangian/intermediate/submodels/Kinematic/CollisionModel/PairCollision/PairModel/PairSpringSliderDashpot/PairSpringSliderDashpot.C

and Line 71-75

Code:

OpenFOAM-2.3.x/src/lagrangian/intermediate/submodels/Kinematic/CollisionModel/PairCollision/PairModel/PairSpringSliderDashpot/PairSpringSliderDashpot.H

we can choose to have adhesive or non-adhesive flows by setting cohesionEnergyDensity greater than 0 or 0 respectively. The surface energy density is a concept used in JKR adhesive contact models. I think by including a value greater than 0 in the code below, DPMFoam will take into account the van der waals forces between contacts. There's a formulae where one can calculate the surface energy as a function of minimum seperation distance of two particle surfaces in the contact region. I will now try and alter this value and analyse the simulation results.

Code:

pairSpringSliderDashpotCoeffs
        {
            useEquivalentSize   no;
            alpha               0.02;
            b                   1.5;
            mu                  0.10;
            cohesionEnergyDensity 0;
            collisionResolutionSteps 12;
        };

But i'm not sure if that cohesion model is based on JKR or an equivalent (linear) cohesion model. I'll investigate this case further

..Not sure if you've tried adjusting this parameter in your numerical simulations?

Thank you very much for your time!
Kevin

Ameya T · July 8, 2015, 04:22

Hi Benz25,

Even I am a new Open foam user. I would like to know that in your example you have mentioned number of particles to be 500, yet in the kinematicCloudPositions file, you have given only one value of position i.e
(0.0005 0.0018 0)
I want to know how will open foam assign co-ordinates to remaining 499 particles?

benz25 · July 18, 2015, 06:35

Hi Ameya T,

I deleted the KinematicCloudPositions file for my simulations and also instead of using manualInjection, I used patchInjection (as ansubru has mentioned in this thread). I defined the area of injection at the inlet using patchInjection. The settings for the paricle injections is shown below. I used patch injection to generate particles at random locations along the walls of the inlet. SOI is the start of injection time during simulation. A linear cohesion model (cohesion energy density) is the only model available in DPM/MPPICFoam to account for particle adhesiveness that is predominant in submicron particles. This could be deployed for simulation. Anyways, hope this helps. Any queries let me know.

Benz

Code:

injectionModels
    {
        model1
        {
            type            patchInjection;
            massTotal       0;
            SOI             0.2;
            parcelBasisType fixed;
            patchName       inlet;
            duration        1;
            parcelsPerSecond  50;
	    nParticle 1;
            U0              (0.5 0 0);
            flowRateProfile constant 1;
            sizeDistribution
            {
                type        fixedValue;
                fixedValueDistribution
                {
                    value   0.001;
                }
            }
        }
    }

Ameya T · August 21, 2015, 16:44

Thanks Kevin for the reply. One more thing which I want to know is would you still be able to use DPMfoam for solid particles in liquid phase?
Also how would you specify number of particles for patchInjection? As in your example you have just kept it as 1.

Ameya T · August 22, 2015, 02:54

Also I would like to know why fixedFluxPressure BCs are applied? Normally I would apply zeroGradient to walls and inlet and fixed value at the outlet. Would this not work?

benz25 · August 31, 2015, 06:43

Hello Ameya,

Unlike the p file in the simpleFOAM solver, the p values for DPMFoam have the phi term embedded...technically I think it should yield the same results but according to this post it's known to have a better convergence? I went through the post before as I had the same question

See for example posts #5, #6:

HTML Code:

http://www.cfd-online.com/Forums/openfoam-solving/82581-i-need-explanations-about-fixedfluxpressure.html

Basically the Boundary conditions are the same I think, either one should yield same numerical results (at least for the case of laminar flows in simple channels)

Hope that helps

Cheers

Kevin

Ameya T · September 2, 2015, 10:34

Hi Kevin!
Thanks again for the reply.

msman · January 19, 2016, 10:12

I am quite new to OpenFOAM and using DPMFoam. The continuous phase field must essentially be air as fluid or it could be water?

ansubru · January 19, 2016, 10:33

Hi there!!

The continuous phase can be any fluid. The relevant properties are entered in the kinematiccloudproperties file. The forces incorporated would also change, and these will have to be individually turned on using the kinematiccloudproperties file again. For instance, when you have water, history forces become important.

Regards
ansubru

msman · January 19, 2016, 11:06

Thanks for your help kind help.

The properties of the closed rectangular computational domain and particles are: 24mm by 77 mm. Inlet velocity at 10 cm/s and have two outlets. Fluid: water Laminar Flow. Solid Particles are made of silica having diameter of 100nm. For me 4-way coupling is very important.

Is it a good decision to use DPMFoam solver?

ansubru · January 20, 2016, 05:38

Hi again!!

From the very brief description of your system it does look plausible that CFD-DEM codes could work. However, as the size of your particles are smaller (nm order of magnitude), I would advise you to calculate the Knudsen numbers to see if rarefied flow dynamics prevail. If that is the case, then you will not be able to treat the fluid as a continuum, as the particles directly interact with the water molecules, so treatment of the system would be markedly different.

Either ways for four-way coupled systems DPMFoam is a good alternative. But be sure that in your case you have a continuum effect as DPMFoam resolves the surrounding fluid as a continuum.

Kind regards

ansubru

msman · January 27, 2016, 10:19

Hello, need help again!!!

I have calculated the Knudsen number and its value is less than 0.01m, So continuum approach is valid. In addition i have made modifications according to my case but now I am getting the following error,

usman@dionysius:~/OpenFOAM/usman-2.3.x/run/lagrangian_src/ZigZagSpacer_Test> #0 Foam::error:rintStack(Foam::Ostream&) in "/nfs/opt/OpenFOAM/OpenFOAM-2.3.x/platforms/linux64GccDPOpt/lib/libOpenFOAM.so"
#1 Foam::sigFpe::sigHandler(int) in "/nfs/opt/OpenFOAM/OpenFOAM-2.3.x/platforms/linux64GccDPOpt/lib/libOpenFOAM.so"
#2 in "/lib64/libc.so.6"
#3 Foam::symGaussSeidelSmoother::smooth(Foam::word const&, Foam::Field<double>&, Foam::lduMatrix const&, Foam::Field<double> const&, Foam::FieldField<Foam::Field, double> const&, Foam::UPtrList<Foam::lduInterfaceField const> const&, unsigned char, int) in "/nfs/opt/OpenFOAM/OpenFOAM-2.3.x/platforms/linux64GccDPOpt/lib/libOpenFOAM.so"
#4 Foam::symGaussSeidelSmoother::smooth(Foam::Field<d ouble>&, Foam::Field<double> const&, unsigned char, int) const in "/nfs/opt/OpenFOAM/OpenFOAM-2.3.x/platforms/linux64GccDPOpt/lib/libOpenFOAM.so"
#5 Foam::smoothSolver::solve(Foam::Field<double>&, Foam::Field<double> const&, unsigned char) const in "/nfs/opt/OpenFOAM/OpenFOAM-2.3.x/platforms/linux64GccDPOpt/lib/libOpenFOAM.so"
#6 Foam::fvMatrix<Foam::Vector<double> >::solveSegregated(Foam::dictionary const&) in "/nfs/opt/OpenFOAM/OpenFOAM-2.3.x/platforms/linux64GccDPOpt/bin/DPMFoam"
#7 Foam::fvMatrix<Foam::Vector<double> >::solve(Foam::dictionary const&) in "/nfs/opt/OpenFOAM/OpenFOAM-2.3.x/platforms/linux64GccDPOpt/bin/DPMFoam"
#8 Foam::fvMatrix<Foam::Vector<double> >::solve() in "/nfs/opt/OpenFOAM/OpenFOAM-2.3.x/platforms/linux64GccDPOpt/bin/DPMFoam"
#9
in "/nfs/opt/OpenFOAM/OpenFOAM-2.3.x/platforms/linux64GccDPOpt/bin/DPMFoam"
#10 __libc_start_main in "/lib64/libc.so.6"
#11
at /home/abuild/rpmbuild/BUILD/glibc-2.14.1/csu/../sysdeps/x86_64/elf/start.S:116

I couldn't figure out which kind of error is it.

ansubru · January 28, 2016, 04:30

Hi again!!

Well, it very hard to try and catch the error you have obtained using the stack you have provided. Some general tips maybe -
1) From the first error in your stack, we could maybe speculate that you are trying to access an unset object. So recheck your setup - geometry, meshing and problem setup including bc's. I would pay particular attention to the geometry and mesh.
2) If you need more help, you may need to furnish some more details of your case (i.e. what are the changes you have made with respect to a working tutorial setup)
3) If you are running the case for the first time, I would suggest you to run OF in debug mode , in that way you are actually pointed to the file and line no. where the first and subsequent errors occur. This is most useful in fine tuning and better understanding how your case behaves.

I hope that helps.

Happy debugging

ansubru

March 25, 2015, 08:40	kinematiccloudproperties	#6
benz25 New Member Join Date: Feb 2015 Posts: 18 Rep Power: 11	Dear Ananda, Thank you very much for your assistance! I went through the websites you posted and took note of the meshing procedure in OF. By following the right hand rule, and using patchinjection, I was able to execute the solver without any issues. I'm going to try cone injection model and try and modify particle properties using different size distribution model and see how it goes... However, in the KinematicCloudProperties file, I'm still not sure what the following terms mean and I'm not sure how or under what situations should I modify those terms (highlighted in red). DPMFoam uses hard sphere model and Hertzian contact if i'm not mistaken?..According to http://www.cfd-online.com/Forums/ope...e-dpmfoam.html does this mean a non adhesive Hertzian contact model is used or adhesive contact i.e. JKR/DMT? Can we select other contact models? And what is the difference between using the openfoam solver DPMFoam and CFDEM (liggghts+openfoam) when it comes to simulating particle-fluid flows? According to website, both codes have 4-way coupling functionality i.e. particle-particle and particle-fluid interactions. So i'm not sure if one of them has something the other doesn't have? Any distinct difference between them? Many thanks again for your assistance!!! Kevin Code: pairCollisionCoeffs { maxInteractionDistance 0.000005; writeReferredParticleCloud no; pairModel pairSpringSliderDashpot; pairSpringSliderDashpotCoeffs { useEquivalentSize no; alpha 0.02; b 1.5; mu 0.10; cohesionEnergyDensity 0; collisionResolutionSteps 12; }; wallModel wallSpringSliderDashpot; wallSpringSliderDashpotCoeffs { useEquivalentSize no; collisionResolutionSteps 12; youngsModulus 1e4; poissonsRatio 0.23; alpha 0.01; b 1.5; mu 0.09; cohesionEnergyDensity 0; }; UName U.air; } Last edited by benz25; March 25, 2015 at 08:53. Reason: contact models

July 18, 2015, 06:35	Particle Injection	#10
benz25 New Member Join Date: Feb 2015 Posts: 18 Rep Power: 11	Hi Ameya T, I deleted the KinematicCloudPositions file for my simulations and also instead of using manualInjection, I used patchInjection (as ansubru has mentioned in this thread). I defined the area of injection at the inlet using patchInjection. The settings for the paricle injections is shown below. I used patch injection to generate particles at random locations along the walls of the inlet. SOI is the start of injection time during simulation. A linear cohesion model (cohesion energy density) is the only model available in DPM/MPPICFoam to account for particle adhesiveness that is predominant in submicron particles. This could be deployed for simulation. Anyways, hope this helps. Any queries let me know. Benz Code: injectionModels { model1 { type patchInjection; massTotal 0; SOI 0.2; parcelBasisType fixed; patchName inlet; duration 1; parcelsPerSecond 50; nParticle 1; U0 (0.5 0 0); flowRateProfile constant 1; sizeDistribution { type fixedValue; fixedValueDistribution { value 0.001; } } } } Farid likes this.

August 31, 2015, 06:43	fixedfluxpressure	#13
benz25 New Member Join Date: Feb 2015 Posts: 18 Rep Power: 11	Hello Ameya, Unlike the p file in the simpleFOAM solver, the p values for DPMFoam have the phi term embedded...technically I think it should yield the same results but according to this post it's known to have a better convergence? I went through the post before as I had the same question See for example posts #5, #6: HTML Code: http://www.cfd-online.com/Forums/openfoam-solving/82581-i-need-explanations-about-fixedfluxpressure.html Basically the Boundary conditions are the same I think, either one should yield same numerical results (at least for the case of laminar flows in simple channels) Hope that helps Cheers Kevin

January 19, 2016, 10:12	Hi,	#15
msman New Member Muhammad Usman Join Date: Nov 2015 Location: Germay Posts: 15 Rep Power: 11	I am quite new to OpenFOAM and using DPMFoam. The continuous phase field must essentially be air as fluid or it could be water?

Similar Threads
Thread	Thread Starter	Forum	Replies	Last Post
Inviscid Drag at subsonic, subcritical Mach #	Axel Rohde	Main CFD Forum	1	November 19, 2001 13:19
SIMPLE method for 2-D flow	N.Natiq	Main CFD Forum	0	June 2, 2001 09:26
flow over a 2D cyl using SIMPLE	T	Main CFD Forum	1	January 27, 2001 08:32
fluid flow fundas	ram	Main CFD Forum	5	June 17, 2000 22:31
momentum under-relaxation for compressible flow with SIMPLE	Mihai ARGHIR	Main CFD Forum	0	April 7, 2000 05:58

March 16, 2015, 05:10		#3
ansubru Member Ananda Kannan Join Date: Feb 2014 Location: Göteborg, Sweden Posts: 55 Rep Power: 12	Hi!! From what i can see, you dont have any problems with particle ínjection (if you did, your simulation will crash at the beginning)... The reason for your crash is that some of the bc's are incorrectly defined (especially the p bc's) as indicated by the error flags.. Code: From function solution::solverDict(const word&) in file matrices/solution/solution.C at line 369 Lookup solver for U.airFinal #0 Foam::error::printStack(Foam::Ostream&) at ??:? #1 Foam::sigFpe::sigHandler(int) at ??:? #2 in "/lib64/libc.so.6" #3 Foam::symGaussSeidelSmoother::smooth(Foam::word const&, Foam::Field<double>&, Foam::lduMatrix const&, Foam::Field<double> const&, Foam::FieldField<Foam::Field, double> const&, Foam::UPtrList<Foam::lduInterfaceField const> const&, unsigned char, int) at ??:? #4 Foam::symGaussSeidelSmoother::smooth(Foam::Field<double>&, Foam::Field<double> const&, unsigned char, int) const at ??:? #5 Foam::smoothSolver::solve(Foam::Field<double>&, Foam::Field<double> const&, unsigned char) const at ??:? #6 Foam::fvMatrix<Foam::Vector<double> >::solveSegregated(Foam::dictionary const&) at ??:? #7 Foam::fvMatrix<Foam::Vector<double> >::solve(Foam::dictionary const&) at ??:? #8 Foam::SolverPerformance<double> Foam::solve<Foam::Vector<double> >(Foam::tmp<Foam::fvMatrix<Foam::Vector<double> > > const&) at ??:? #9 at ??:? You cannot define a 'pressureinletoutlet' or ínterstitialinletoutlet' boundary as a fixedfluxpressure boundary... in such case rename the corresponding p boundary as 'fixedValue'... this might solve your problem... Additionally, i can see that your particle size is 0.00005m, hence correspondingly the maxinteractiondistance defined in the kinematiccloudproperties file must also be 0.00005. Additionally, if you want better control over your simulations, i would advice you to setup your case from scratch,i.e geometry and meshing, in OF (your setup is relatively simple)... In this way you are sure that the solutions are not affected by the different ways in which your case has been meshed... Hope this msg will help you Kind regards ansubru

March 19, 2015, 09:37		#5
ansubru Member Ananda Kannan Join Date: Feb 2014 Location: Göteborg, Sweden Posts: 55 Rep Power: 12	Hi Kevin!! I have a few suggestions for you : --> Firstly run the command çheckMesh' to see mesh statistics... I have a feeling you have not followed the 'right hand rule' while assigning patches and that your mesh is quite poor (please read the blockMesh help in detail before you begin meshing)... refer here..http://www.openfoam.org/docs/user/blockMesh.php --> Secondly, you need to define your injection correctly, if your mesh has been 'drawn' in a different co-ordinate system, you could have issues when injections have to be initialized... And try to use alternate particle injection strategies.. I use patchinjection, and that seems to work quite well for me ... --> Thirdly, try not to use the exact young's modulus while computing... Use a value 1000 times lower (for example try not to exceed the value used in the DPMFoam tutorial)... Your simulation will still be correct, dont worry.. Its an approximation almost all DEM works use so as to ease the computational burden... I have commented on a post regarding the same.. read more here.. http://www.cfd-online.com/Forums/ope...tml#post513475 Hope that helped!! Regards ansubru

March 25, 2015, 11:19		#7
ansubru Member Ananda Kannan Join Date: Feb 2014 Location: Göteborg, Sweden Posts: 55 Rep Power: 12	Hi Kevin!! Firstly, DPM-Foam is a 'soft sphere' contact model based solver. The entire concept of CFD-DEM has been developed by Tsuji et.al (look at the paper below) while studying plug flow conveying. All relevant parameters for the method of DEM are detailed in their paper, you can have a detailed read and I am sure that will answer all your questions... https://www.dropbox.com/sh/6bis355hp...ihf45kyra?dl=0 Equivalent size is a concept that comes into play when you have non spherical particles. My suggestion would be for you to read about hertzian contact models and why they are universally accepted as one of the best contact models available for simulating multi-particle contact. Additionally, the DPMFoam solver can also account for cohension, but I am not very aware of the models behind that. The config file you have shared is for 'non-adhesive' contacts. As of now, DPMFoam only uses linear and hertzian contact pair spring slider dashpot models for contact. Any new model, will unfortunately have to be developed you... And finally, about the differences between the various solvers (i.e CFDEM-LIGGGHTS, etc.) there is'nt much. They have the same features and are equally effective in resolving granular flow physics. LIGGGHTS has one major advantage that there are lot of supporting documentation/tutorials on setting up and analyzing a simulation. Apart from that, there is also the added capability of modelling dynamic meshes- moving geometries (which unfortunately cant be done with DPMFoam right now)... However, DPMFoam is more malleable to your needs than LIGGGHTS, i.e if you of course read and understand the source code. Else both tools are equally capable and its a matter of choice really. Personally, I couldnt get LIGGGHTS to work some incompatibility with my hardware, hence I adopted DPMFoam Hope that helped All the best in your quest!! ansubru rmn_990 likes this.

April 2, 2015, 12:01		#8
benz25 New Member Join Date: Feb 2015 Posts: 18 Rep Power: 11	Hi Ananda, Thank you very much for your reply. I did have some compatibility issues with liggghts as well. After spending some time with the installation, I did manage to install it using Ubuntu 12 LTS and OpenFOAM 2.1.1. I used either Helyx or Salome for meshing. Thanks for sharing the document pertaining to CFD-DEM developed by Tsuji (1993). I've read two of Tsuji's papers. I see they use the spring dashpot model for the forces whilst the Hertz theory is used to evaluate the stiffness. So DPMFoam allows users to select either linear-spring–dashpot model or non-linear Hertzian-spring–dashpot by setting b=1 or 1.5 respectively: Code: [OpenFOAM-2.3.x/src/lagrangian/intermediate/submodels/Kinematic/CollisionModel/PairCollision/PairModel/PairSpringSliderDashpot/PairSpringSliderDashpot.H Interestingly, after going through additional papers regarding surface energies (used in cohesion models), and then looking back at the DPMFoam code, I noted that in Line 217 Code: OpenFOAM-2.3.x/src/lagrangian/intermediate/submodels/Kinematic/CollisionModel/PairCollision/PairModel/PairSpringSliderDashpot/PairSpringSliderDashpot.C and Line 71-75 Code: OpenFOAM-2.3.x/src/lagrangian/intermediate/submodels/Kinematic/CollisionModel/PairCollision/PairModel/PairSpringSliderDashpot/PairSpringSliderDashpot.H we can choose to have adhesive or non-adhesive flows by setting cohesionEnergyDensity greater than 0 or 0 respectively. The surface energy density is a concept used in JKR adhesive contact models. I think by including a value greater than 0 in the code below, DPMFoam will take into account the van der waals forces between contacts. There's a formulae where one can calculate the surface energy as a function of minimum seperation distance of two particle surfaces in the contact region. I will now try and alter this value and analyse the simulation results. Code: pairSpringSliderDashpotCoeffs { useEquivalentSize no; alpha 0.02; b 1.5; mu 0.10; cohesionEnergyDensity 0; collisionResolutionSteps 12; }; But i'm not sure if that cohesion model is based on JKR or an equivalent (linear) cohesion model. I'll investigate this case further..Not sure if you've tried adjusting this parameter in your numerical simulations? Thank you very much for your time! Kevin rmn_990 likes this.

July 8, 2015, 04:22		#9
Ameya T New Member Join Date: Jul 2015 Posts: 17 Rep Power: 11	Hi Benz25, Even I am a new Open foam user. I would like to know that in your example you have mentioned number of particles to be 500, yet in the kinematicCloudPositions file, you have given only one value of position i.e (0.0005 0.0018 0) I want to know how will open foam assign co-ordinates to remaining 499 particles?

August 21, 2015, 16:44		#11
Ameya T New Member Join Date: Jul 2015 Posts: 17 Rep Power: 11	Thanks Kevin for the reply. One more thing which I want to know is would you still be able to use DPMfoam for solid particles in liquid phase? Also how would you specify number of particles for patchInjection? As in your example you have just kept it as 1.

August 22, 2015, 02:54		#12
Ameya T New Member Join Date: Jul 2015 Posts: 17 Rep Power: 11	Also I would like to know why fixedFluxPressure BCs are applied? Normally I would apply zeroGradient to walls and inlet and fixed value at the outlet. Would this not work?

September 2, 2015, 10:34		#14
Ameya T New Member Join Date: Jul 2015 Posts: 17 Rep Power: 11	Hi Kevin! Thanks again for the reply.

January 19, 2016, 10:33		#16
ansubru Member Ananda Kannan Join Date: Feb 2014 Location: Göteborg, Sweden Posts: 55 Rep Power: 12	Hi there!! The continuous phase can be any fluid. The relevant properties are entered in the kinematiccloudproperties file. The forces incorporated would also change, and these will have to be individually turned on using the kinematiccloudproperties file again. For instance, when you have water, history forces become important. Regards ansubru

January 19, 2016, 11:06		#17
msman New Member Muhammad Usman Join Date: Nov 2015 Location: Germay Posts: 15 Rep Power: 11	Thanks for your help kind help. The properties of the closed rectangular computational domain and particles are: 24mm by 77 mm. Inlet velocity at 10 cm/s and have two outlets. Fluid: water Laminar Flow. Solid Particles are made of silica having diameter of 100nm. For me 4-way coupling is very important. Is it a good decision to use DPMFoam solver?

January 20, 2016, 05:38		#18
ansubru Member Ananda Kannan Join Date: Feb 2014 Location: Göteborg, Sweden Posts: 55 Rep Power: 12	Hi again!! From the very brief description of your system it does look plausible that CFD-DEM codes could work. However, as the size of your particles are smaller (nm order of magnitude), I would advise you to calculate the Knudsen numbers to see if rarefied flow dynamics prevail. If that is the case, then you will not be able to treat the fluid as a continuum, as the particles directly interact with the water molecules, so treatment of the system would be markedly different. Either ways for four-way coupled systems DPMFoam is a good alternative. But be sure that in your case you have a continuum effect as DPMFoam resolves the surrounding fluid as a continuum. Kind regards ansubru

January 27, 2016, 10:19		#19
msman New Member Muhammad Usman Join Date: Nov 2015 Location: Germay Posts: 15 Rep Power: 11	Hello, need help again!!! I have calculated the Knudsen number and its value is less than 0.01m, So continuum approach is valid. In addition i have made modifications according to my case but now I am getting the following error, usman@dionysius:~/OpenFOAM/usman-2.3.x/run/lagrangian_src/ZigZagSpacer_Test> #0 Foam::error:rintStack(Foam::Ostream&) in "/nfs/opt/OpenFOAM/OpenFOAM-2.3.x/platforms/linux64GccDPOpt/lib/libOpenFOAM.so" #1 Foam::sigFpe::sigHandler(int) in "/nfs/opt/OpenFOAM/OpenFOAM-2.3.x/platforms/linux64GccDPOpt/lib/libOpenFOAM.so" #2 in "/lib64/libc.so.6" #3 Foam::symGaussSeidelSmoother::smooth(Foam::word const&, Foam::Field<double>&, Foam::lduMatrix const&, Foam::Field<double> const&, Foam::FieldField<Foam::Field, double> const&, Foam::UPtrList<Foam::lduInterfaceField const> const&, unsigned char, int) in "/nfs/opt/OpenFOAM/OpenFOAM-2.3.x/platforms/linux64GccDPOpt/lib/libOpenFOAM.so" #4 Foam::symGaussSeidelSmoother::smooth(Foam::Field<d ouble>&, Foam::Field<double> const&, unsigned char, int) const in "/nfs/opt/OpenFOAM/OpenFOAM-2.3.x/platforms/linux64GccDPOpt/lib/libOpenFOAM.so" #5 Foam::smoothSolver::solve(Foam::Field<double>&, Foam::Field<double> const&, unsigned char) const in "/nfs/opt/OpenFOAM/OpenFOAM-2.3.x/platforms/linux64GccDPOpt/lib/libOpenFOAM.so" #6 Foam::fvMatrix<Foam::Vector<double> >::solveSegregated(Foam::dictionary const&) in "/nfs/opt/OpenFOAM/OpenFOAM-2.3.x/platforms/linux64GccDPOpt/bin/DPMFoam" #7 Foam::fvMatrix<Foam::Vector<double> >::solve(Foam::dictionary const&) in "/nfs/opt/OpenFOAM/OpenFOAM-2.3.x/platforms/linux64GccDPOpt/bin/DPMFoam" #8 Foam::fvMatrix<Foam::Vector<double> >::solve() in "/nfs/opt/OpenFOAM/OpenFOAM-2.3.x/platforms/linux64GccDPOpt/bin/DPMFoam" #9 in "/nfs/opt/OpenFOAM/OpenFOAM-2.3.x/platforms/linux64GccDPOpt/bin/DPMFoam" #10 __libc_start_main in "/lib64/libc.so.6" #11 at /home/abuild/rpmbuild/BUILD/glibc-2.14.1/csu/../sysdeps/x86_64/elf/start.S:116 I couldn't figure out which kind of error is it.

January 28, 2016, 04:30		#20
ansubru Member Ananda Kannan Join Date: Feb 2014 Location: Göteborg, Sweden Posts: 55 Rep Power: 12	Hi again!! Well, it very hard to try and catch the error you have obtained using the stack you have provided. Some general tips maybe - 1) From the first error in your stack, we could maybe speculate that you are trying to access an unset object. So recheck your setup - geometry, meshing and problem setup including bc's. I would pay particular attention to the geometry and mesh. 2) If you need more help, you may need to furnish some more details of your case (i.e. what are the changes you have made with respect to a working tutorial setup) 3) If you are running the case for the first time, I would suggest you to run OF in debug mode , in that way you are actually pointed to the file and line no. where the first and subsequent errors occur. This is most useful in fine tuning and better understanding how your case behaves. I hope that helps. Happy debugging ansubru