[dussa]

Hi everyone,

I am having a bit of trouble with implementing my own geometry and case in DPMFoam. I have taken the DPMFoam/Goldschmidt tutorial case and modified it so that I can investigate turbulent flow of water through a buoyant bed of particles.

The things that I have made changes to (most of which probably have no consequence with the error I am having issues with) are:

• the geometry with fluent3DMeshToFoam

• the continuous phase from air to water

• patch definitions for walls and inlet/outlet in 0/*

• different particle injection file consistent with new geometry

I have initially run the simulation keeping it as laminar to ensure that some things such as buoyancy of the particles and inlet/outlets are set correctly, which they seem to be. I have now tried changing from laminar to k-epsilon RAS turbulence model, which is where my issues seem to be starting.

I have changed the continuousPhaseName air ---> continuousPhaseName water in the constant/transportProperties while also adding the required files in the /0 directory and altering some names due to the use of water in my naming such that my case directory has the following structure now:

TMSR-SF1
├── 0
│** ├── p
│** ├── U.water
│** ├── epsilon.water
│** ├── k.water
│** └── nut.water
├── constant
│** ├── g
│** ├── kinematicCloudPositions
│** ├── kinematicCloudProperties
│** ├── polyMesh/...
│** ├── transportProperties
│** └── turbulenceProperties.water
├── system
│** ├── controlDict
│** ├── fvSchemes
│** └── fvSolution
└── TMSRSF1ascii.msh

I was getting various errors when trying to run the case, and having done some research into them found that it was mainly to with naming problems in the other case files (leaving .water off many of them). That got me as far as this now. I am now getting the following error:

/*---------------------------------------------------------------------------*\
| =========                 |                                                 |
| \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox           |
|  \\    /   O peration     | Version:  4.x                                   |
|   \\  /    A nd           | Web:      www.OpenFOAM.org                      |
|    \\/     M anipulation  |                                                 |
\*---------------------------------------------------------------------------*/
Build  : 4.x-b5331fdc8656
Exec   : DPMFoam
Date   : Dec 19 2017
Time   : 12:22:57
Host   : "XPS15"
PID    : 15604
Case   : /home/robert/OpenFOAM/robert-4.x/run/TMSR-SF1-turbulent
nProcs : 1
sigFpe : Enabling floating point exception trapping (FOAM_SIGFPE).
fileModificationChecking : Monitoring run-time modified files using timeStampMaster
allowSystemOperations : Allowing user-supplied system call operations

// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
Create time

Create mesh for time = 0


PIMPLE: Operating solver in PISO mode


Reading g
Reading field U

Reading field p

Reading/calculating continuous-phase face flux field phic

Creating turbulence model

Selecting incompressible transport model Newtonian
Creating field alphac

Constructing kinematicCloud kinematicCloud
Constructing particle forces
    Selecting particle force ErgunWenYuDrag
    Selecting particle force gravity
    Selecting particle force pressureGradient
Constructing cloud functions
    none
Constructing particle injection models
Creating injector: model1
Selecting injection model manualInjection
    Constructing 3-D injection
    Choosing nParticle to be a fixed value, massTotal variable now does not determine anything.
Selecting distribution model fixedValue
Selecting dispersion model none
Selecting patch interaction model localInteraction
    Interaction fields will not be written
Selecting stochastic collision model none
Selecting surface film model none
Selecting U integration scheme Euler
Selecting collision model pairCollision
Selecting pair model pairSpringSliderDashpot
Selecting wall model wallSpringSliderDashpot
Building InteractionLists with interaction distance 0.0025
    Building referred interaction lists
    Building direct interaction lists
Selecting turbulence model type RAS
Selecting RAS turbulence model kEpsilon
kEpsilonCoeffs
{
    Cmu             0.09;
    C1              1.44;
    C2              1.92;
    C3              0;
    sigmak          1;
    sigmaEps        1.3;
}


Starting time loop

Courant Number mean: 1.67985e-07 max: 9.86704e-05
Time = 2e-05

Evolving kinematicCloud

Solving 3-D cloud kinematicCloud

Cloud: kinematicCloud injector: model1
    Added 100 new parcels

Cloud: kinematicCloud
    Current number of parcels       = 100
    Current mass in system          = 0.706858
    Linear momentum                 = (0 0 0.000277371)
   |Linear momentum|                = 0.000277371
    Linear kinetic energy           = 5.44201e-08
    model1:
        number of parcels added     = 100
        mass introduced             = 0.706858
    Parcel fate (number, mass)      : patch inlet
      - escape                      = 0, 0
      - stick                       = 0, 0
    Parcel fate (number, mass)      : patch walls
      - escape                      = 0, 0
      - stick                       = 0, 0
    Parcel fate (number, mass)      : patch outlet
      - escape                      = 0, 0
      - stick                       = 0, 0
    Rotational kinetic energy       = 0

PIMPLE: iteration 1
smoothSolver:  Solving for U.waterx, Initial residual = 1, Final residual = 3.58584e-15, No Iterations 1
smoothSolver:  Solving for U.watery, Initial residual = 1, Final residual = 3.97685e-15, No Iterations 1
smoothSolver:  Solving for U.waterz, Initial residual = 1, Final residual = 1.15756e-16, No Iterations 1
GAMG:  Solving for p, Initial residual = 1, Final residual = 0.00991846, No Iterations 13
time step continuity errors : sum local = 3.36373e-09, global = 4.33885e-10, cumulative = 4.33885e-10
GAMG:  Solving for p, Initial residual = 0.605884, Final residual = 8.04761e-07, No Iterations 17
time step continuity errors : sum local = 3.28417e-12, global = -4.22639e-13, cumulative = 4.33462e-10
No finite volume options present



--> FOAM FATAL ERROR: 

    request for surfaceScalarField phi from objectRegistry region0 failed
    available objects of type surfaceScalarField are

6
(
alphaPhic
(interpolate((((1|A(U.water))*cloudVolSUSu)|rho.water))+(Dp*(g&S)))
phi.water_0
phi.water
Dp
alphacf
)


    From function const Type& Foam::objectRegistry::lookupObject(const Foam::word&) const [with Type = Foam::GeometricField<double, Foam::fvsPatchField, Foam::surfaceMesh>]
    in file /home/robert/OpenFOAM/OpenFOAM-4.x/src/OpenFOAM/lnInclude/objectRegistryTemplates.C at line 193.

FOAM aborting

#0  Foam::error::printStack(Foam::Ostream&) at ??:?
#1  Foam::error::abort() at ??:?
#2  Foam::GeometricField<double, Foam::fvsPatchField, Foam::surfaceMesh> const& Foam::objectRegistry::lookupObject<Foam::GeometricField<double, Foam::fvsPatchField, Foam::surfaceMesh> >(Foam::word const&) const at ??:?
#3  Foam::turbulentMixingLengthDissipationRateInletFvPatchScalarField::updateCoeffs() at ??:?
#4  Foam::GeometricField<double, Foam::fvPatchField, Foam::volMesh>::Boundary::updateCoeffs() at ??:?
#5  Foam::RASModels::kEpsilon<Foam::PhaseIncompressibleTurbulenceModel<Foam::singlePhaseTransportModel> >::correct() at ??:?
#6  ? at ??:?
#7  __libc_start_main in "/lib/x86_64-linux-gnu/libc.so.6"
#8  ? at ??:?
Aborted (core dumped)

It creates the fields, inserts the particles, and begins running the first PIMPLE iteration and then spits out the error:

No finite volume options present



--> FOAM FATAL ERROR: 

    request for surfaceScalarField phi from objectRegistry region0 failed
    available objects of type surfaceScalarField are

6
(
alphaPhic
(interpolate((((1|A(U.water))*cloudVolSUSu)|rho.water))+(Dp*(g&S)))
phi.water_0
phi.water
Dp
alphacf
)


    From function const Type& Foam::objectRegistry::lookupObject(const Foam::word&) const [with Type = Foam::GeometricField<double, Foam::fvsPatchField, Foam::surfaceMesh>]
    in file /home/robert/OpenFOAM/OpenFOAM-4.x/src/OpenFOAM/lnInclude/objectRegistryTemplates.C at line 193.

FOAM aborting

I have read the OpenFOAM Wiki on what phi is, but I still am not sure. I have tried changing all of the instances of phi to phi.water (as I did with other things such as U U.water etc.) but nothing changes. I also can't find where there are any surfaceScalarFields defined in this problem.

Is there someone that can point me in the correct direction? I have tried writing my own 0/phi.water file but that also didn't seem to help (mainly because I was not sure what information is actually required in the phi.water file).

Please let me know if anyone has any suggestions, and please find attached the case files.
Kind regards,
Robert

[blebon]

Since DPM and MPPICFoam are mostly similar, I modified your case to mimic the turbulent cyclone tutorial and got it working with LES (attached)

Not sure about the above error. Maybe you did not define turbulence divSchemes properly. Note that the variable names in the 0 directory are also fixed in the header (k -> k.water for example, you changed only the file names but not the headers). Also, you need to add the divSchemes for each solved variables. You did not add div(alphaPhic,k.water) Gauss limitedLinear 1; for example.

In fvSolution, you need to add .water to all other solved variables as well, for example "(U|k|epsilon|omega).water" instead of "(U.water|k|epsilon|omega)" in your original case.

I am also interested in working with Euler-Lagrange solvers, so your case helped me to get DPM working in OpenFOAM. I hope these changes helps you as well.

Regards,

Bruno

[dussa]

Thank you for your reply Bruno! And thank you for your upload showing the changes you made to allow the simulation to work with LES. I am still puzzled then as to why it won't work with RAS.

I ran your simulation, and although it runs, the motion of the pebbles is not as expected (they seem to hit walls before reaching the top of the domain/some get stuck middle of the domain), perhaps I will have to go back a bit further now and reduce the interaction forces and increase relaxation factors for the fluid-solid phase transfers.

Since making the first post I actually went a step further back and tried to simply simulate a laminar model. I made a post about it in another thread but I'll copy the relevant stuff here.

Using a shortened Goldschmidt geometry with fewer particles and no inlet velocity:

When the simulation is with air (1.2 kg/m3) the particles (rho0 = 2526 kg/m3) fall to the bottom of the geometry as expected and quickly lose most of their kinetic energy forming a layer with very little movement.

When I try this with the case altered by using water (rho = 1000 kg/m3) and reducing the particle density (rho0 = 750 kg/m3) I find that the particles move in a very erratic manner, and do not simply float from their injection location to the top of the fluid domain (you can even see that after 0.01sec, the first image, they have already moved in a strange way and are not in a straight line). They seem to be attracted to the walls/symmetry boundaries, and away from the centre of the fluid domain as you can see in the attached images after 1 second of physical time.

So the issue seems to be with the pebbles being in a more dense phase, not so much the introduction of turbulence, as the simulation with air as the working fluid produces realistic motion of the pebbles. I am wondering if you or anyone else has any ideas on this?

Kind regards,
Robert

[dussa]

Case directory for the shortened Goldschmidt geometry with water as fluid, less dense particles, solver as laminar.

[procEng]

We implemented a geometry into the MPPICFoam-solver and did some modifications with the turbulenceProperties. After this we run into the same error. What we found out: Depending on the boundarys in your case and the turbulence model in use, you have to declare phi,k etc. in multiple Fields (k,epsilon) in the 0 Directory.

Durs

[dussa]

Thanks you for you reply Durs,
My understanding of what you have described is that I need to declare in $CASE_DIR/0 the fields that I see in the appearing in the later time directories?

Looking at my time directories they have the following contents:

goldschmidtTestBuoyant/
├── 0/
│** ├── p
│** ├── U.h2o
├── 0.01/
│** ├── lagrangian/
│** ├── uniform/
│** ├── p
│** ├── U.h2o
│** ├── alpha.h2o
│** ├── mu.h2o
│** ├── phi.h2o
│** ├── rho.h2o
│** ├── kinematicCloud:UCoeff
│** └── kinematicCloud:UTrans
.
.
.
etc.

I tried looking in the new files that are created in the 0.01 and onward directories to get an idea as to how to implement them in /0 but they are (half) unreadable binary (EDIT: I have since changed writeFormat from binary to ascii in the controlDict), so I am not entirely sure how to initially implement them in /0.

My simulation currently is laminar, so I have no need for k, epsilon, etc.

I will have a go and report back if I make any progress.

Regards,
Rob

[procEng]

Hey Rob

hm, when i read my post, its pretty unclear to my self what i have tried to say:

If the error message shows phi: it can help to check the boundaryFields in all the files of the 0 Directory. In our case there are two Inlets and a Outlet and after we added the red lines:

 ... boundaryField
inlet1
{
type ...
ref Value...
phi    phi.air;
}

inlet2
{
type ...
ref Value...
phi    phi.air;
}

outlet
{
type ...
inletValue...
value....
phi    phi.air;
}
...

to this three boundarys (wall didn't need the line) the error message changed from:

request for volScalarField phi....

to

request for volScalarField k...

so, we did the same for k and the errors were gone.

[dussa]

That seemed to work! Thank you very much for pointing out that I was missing some declarations in the volScalarFields. I also uncovered a couple other issues whilst I was working on fixing this up, namely that in fvSolution I had omega instead of epsilon (I guess the case I was adapting was k-omega and I never noticed, mainly because there was never an /0/omega file). So I had to just change:

"(U.water|k.water|omega.water)"
.....
"(U.water|k.water|omega.water)Final"

to

"(U.water|k.water|epsilon.water)"
.....
"(U.water|k.water|epsilon.water)Final"

I also went through and found that I had in some cases inserted the phi phi.water in too many places in the volVectorField and volScalarField in /0, so I commented them out. For example, in /0/U.water I found that the phi term is only needed at the outlet patch, and not the inlet:

FoamFile
{
    version     2.0;
    format      binary;
    class       volVectorField;
    location    "0";
    object      U.water;
}
// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
dimensions      [0 1 -1 0 0 0 0];
internalField   uniform (0 0 0);
boundaryField
{
    walls
    {
        type            noSlip;
    }
    outlet
    {
        type            pressureInletOutletVelocity;
        phi             phi.water;
        inletValue      uniform (0 0 0);
        value           uniform (0 0 0);
    }
    inlet
    {
        type            interstitialInletVelocity;
        inletVelocity   uniform (0 0 1.875);
        value           uniform (0 0 1.875);
//      phi             phi.water; // not required by inlet patch
        alpha           alpha.water;
    }
}

I will run the turbulent simulation. However I am still having issues with just laminar cases where it appears that the momentum transfer from the solid particles to the fluid is too strong, and even in a system with no flow (no inlet or outlet flows) there are huge flow instabilities that create unrealistic flows.

Do you have any experience with how to reduce the severity of the momentum transfer from solid phase to fluid phase? I could just remove the drag model (currently using ErgunWenYu as I have in the past with StarCCM+ simulations) but then I lose the void fraction mapping of the solid phase to the fluid domain, which is what I am interested in.

In my simulation, at its operating state, the pebbles will be floating against the top of the container, with flow passing through the bed from bottom to top, so there will be practically no movement of the pebbles anyway. If there is any method of reducing the effect of the momentum being transferred to the fluid from the pebbles, whilst retaining the void fraction mapping then I would greatly appreciate if someone could shed some light on that.

I will try and have a dig around the source code of the models, but I am relatively green when it comes to programming, and the further I delve into OpenFoams source code, the more and more lost I seem to be getting!

Kind regards,
Rob

[joaran]

Hi Robert, were you finally able to fix the "instabilities" issues? If so, could you please give me a hint?



Regards,

Joa

[dussa]

If I remember correctly my issue was solved by increasing the length of the timestep.


Regards
Rob

[farzadmech]

Hello friends

what I am looking for is to find the Y position of particles(bubbles) and put it in below code;

const KinematicParcel<ParcelType>& p;


//Increase in diameter due to rising of bubble using ideal gas relation P1*V1=P2*V2 ===> V2 = (P1*V1)/P2
scalar V0 =Vinit*(Patm+(Ymax - 0.)*rhoo*gg) / (Patm+(Ymax - p.position()[2])*rhoo*gg);

When I do this, I have got an error which says that the initialization has not done correctly, what should I do?