[mattijs]

Info and Pout buffer output until an 'endl' is written.

AUTO_WRITE should work if the variable still is in scope when runTime.write() gets called.

If not you can always force it to write: e.g. alpaha.write()

[alberto]

I finally understood the problem.
I started to work with a version of the solver in my user directory, then I moved to the version in the OpenFOAM main directory, without removing the libraries generated by the old one.
The linker still used them to link the executable, so this is the reason my changes to the code didn't take effect.
I found this out by commenting out g0, and noticing that gdb still was saying me there was a floating point exception in g0. :-)

Now I'm doing some other test to see if the alternative g0 formulation (also adopted in MFIX) works properly.

As always, thanks for your precious help.

Alberto

[alberto]

Hello, what solver settings were used in the examples you showed me some month ago where you did the comparison between the algebraic and the differential transport equation of the granular energy?

I'm obtaining better results by approximating g0 around the packing limit and using SuperBee and limiting the gradient.

My doubt is what flux do I specify for div((-nuEffa*grad(Ua).T())) and div((-nuEffb*grad(Ub).T())) when using an upwind scheme?

Best regards,
Alberto

[niklas]

Dont remember.

the div-schemes for the neff terms should be linear.

[mattijs]

1) mesh.cellCells()[cellI] gives the neighbours of a cell (see primitiveMesh.H)
2) with those cell labels you can just index into any volScalarField
3) loop over all faces of a cell (mesh.cells()[cellI]) and check if any one is not an internal face (mesh.isInternalFace(faceI))
(see primitiveMesh.H)

or loop over all the boundary faces and mark the owner cell (mesh.faceOwner()[faceI]). Then pick up marked cells. See User Guide about mesh ordering (boundary faces after internal faces)

[olivier]

Hi,
1- I have a problem of convergence with twoEulerFoam. I use the different geometry and I change the schemes and time step. herewith one of my files .msh
I patch the value 1 for alpha in half-top, in Inlet and in outlet in my bubble column. 0 in the half-buttom.
the diameter of Air is 0.003m
the velocity in inlet is:
Ua=0.0212 m/s
Ub=0;
I have a good result with fluent.
with Foam Ur Courant Number grow up and explode.
with or without turbulence.
thank you for your help.

2- how can I attache my case ?

[alberto]

To attach a file, insert the text

\ attach{Description of your file}

without the space after \
The system will ask you to upload the file after you post the message. More information here: http://www.cfd-online.com/cgi-bin/Op...?pg=formatting

Regards,
Alberto

[olivier]

Hi Alberto,
I think the problem comes from the initialization of alpha (volume fraction).

After the resolution of alpha-equation, when alpha<0 I make it =0
and if alpha>1 I make it=1
in this case the divergence is lower.
The results by fluent with the same mesh and same conditions are very good.

I will add the code of population balance witch is already written for fluent (UDF) when I resolve this problem.
thanks a lot for taking time to look to my case.

[olivier]

Hi Alberto,
I think the problem comes from the initialization of alpha (volume fraction).

After the resolution of alpha-equation, when alpha<0 I make it =0
and if alpha>1 I make it=1
in this case the divergence is lower.
The results by fluent with the same mesh and same conditions are very good.

I will add the code of population balance witch is already written for fluent (UDF) when I resolve this problem.
thanks a lot for taking time to look to my case.

[alberto]

Hello Oliver,
I can't download your attached file. Could you email it to me? My e-mail address should be visible in the forum.

However, the limitation of alpha between zero and one is managed in OpenFOAM by rewriting the equation for alpha in a fully conservative form, as proposed by Henry and by using a limited scheme, like limitedLinear01 or similar. Everything is well explained in the PhD thesis of Henry Rusche, downloadable from the site of Hrvoje.

In my calculations, the instability is related to the kinetic theory of granular flow, because of the presence of the radial distribution function which becomes very big if the solids volume fraction gets close to the packing value. This causes instabilities in the momentum equation of the dispersed phase, due to the presence of the granular pressure gradient.

I'm just testing an algorithm to avoid these problems, as I explained in another thread.

Best regards,
Alberto

[frank]

Hi everyone,
I wonder how viscosity include bulk and collisional-kinetic distribution and pressrue of particles of kinetic theory model transfer to N-S equation of solid phase in implement of twoPhaseEulerFoam.
I think it is a bug! I guess this codes were originlly developed for gas-liquids flow using the nature visosity of liquids in N-S equation of liquids.

I managed to revised the source code of UEqns.H as follows,

nuEffa=kineticTheory.mua()/rhoa;
volTensorField Rca = -nuEffa*(fvc::grad(Ua)().T());
Rca += (2.0/3.0)*sqr(Ct)*I*k - kineticTheory.pa()/rhoa*I +kineticTheory.lambda()/rhoa*tr(fvc::grad(Ua)().T())*I - (2.0/3.0)*I*tr(Rca);

Anybody tell me whether my doubt is right.

What are your suggests?
Thanks in advance,

Best regards
Q. Li

[alberto]

Hello,
the granular pressure gradient is included in the pressure equation (pEqn.H):

if (kineticTheory.on())
{
phiDraga -= rUaAf*fvc::snGrad(kineticTheory.pa()/rhoa)*mesh.magSf();
}

You're right about the bulk viscosity, which is not included in the original solver.
To keep the solver general, you can modify UEqn.H by replacing the lines:

volTensorField Rca = -nuEffa*(fvc::grad(Ua)().T());
Rca += (2.0/3.0)*sqr(Ct)*I*k - (2.0/3.0)*I*tr(Rca);

in

volTensorField gradUaT = fvc::grad(Ua)().T();
volTensorField Rca = -nuEffa*gradUaT;
Rca += (2.0/3.0)*sqr(Ct)*I*k - (2.0/3.0)*I*tr(Rca);

if (kineticTheory.on())
{
Rca += kineticTheory.lambda()/rhoa*I*tr(-gradUaT);
}

The granular viscosity is taken into account in twoPhaseEulerFoam.C:

if (kineticTheory.on())
{
kineticTheory.solve();
nuEffa += kineticTheory.mua()/rhoa;
}

Regards,
Alberto

[frank]

Thanks for your exact and minute answer! I have benefited from it.
By the way, I'm also trying to simulate a bubbling bed with a jet, Schematic representation, condition and physical properties of fluidized bed come from D. J. Patil[Critical comparision of hydrodynamic models for gas-solid fluidized beds part I: bubbling gas-solid fluidized beds operated with a jet. Chem. Eng. Sci. 60(2005)57-72]. which I used the drag model of your distribution - GidaspowErgunWenYu, viscosityModel-Gidaspow, conductivtyModel-Gidaspow. But I can't obtain a well-pleasing result. can you give me some advices!
In addition I have commented codes of fictional stress model.

warmly regards,
Q. Li

[alberto]

What kind of problems do you have?
Can you complete the calculation, but the bubble shape is wrong? Or you can't obtain a converged solution?

Regards,
Alberto

[frank]

check your Email or see follows:

just as what you said on Sunday, September 18, 2005 - 08:19's Post
"In my opinion the problems depend on the management of the packing limit.

If I use the particle particle interaction force by setting g0 != 0 in the ppProperties dictionary, I get convergence, but the solid volume fraction is wrong. The bubble has a wrong form and in some case it just breaks instead of closing on itself and going up.

If I set g0 = 0 in order to manage the packing limit with the kinetic theory approach, the result seems to be better, but it's impossible to complete a calculation due to convergence problems. After a few iteration I get a "nan" in the continuity equation residual. This should be related to the radial distribution function g0 which becomes infinite if the packing limit is reached. "

I have encountered the same instance. Are you solve this problem? I hope you can share you experiments with me.

regards
Q. Li

[niklas]

As I undertand this code has gone through some improvements and if you want it included in the official 1.3 version, it might be good to send it to Henry.
Cleaned up of course.

Niklas

[alberto]

Yes. I'm working on cleaning the code right now and on making a list of the changes I did.

I have some test case too which can be of interest, but I haven't extended the FoamX configuration files yet.

Alberto

[alberto]

I've already implemented and generalised that method, which seems to work well, even though it requires a small time step (10^-5).

I'm waiting for a paper acceptance and I'll release it.

In the meanwhile, pay attention that the method presented in my previous post is valid only for cells of equal volume (I extended it to the general case).
Also the method has to be applied not only to the packed cells but also to the neighbouring ones in order not to lose mass.

Regards,
Alberto

[guanghaowu]

I wonder if twoPhaseEulerFoam is suitable for a pipe flow with solid paticles?

I am trying to calculate the following three cases (single processor) with twoPhaseEulerFoam,
and found all of the cases fell in divergence. Due to the gravity force,
the solid particles will be deposited on the bottom of the pipe.
So the alpha gets greater and greater and becomes alpha > 0.6~0.9
( the value is depend on the case) . Then the solver diverges.

â–*Simulation cases
3D pipe flow (liquid) with solid particles

The initial conditions:
alpha = 1e-05 for internalField; 0.2 for inlet boundary
Ua = 0m/sec for internalField; 1.0m/sec for inlet boundary
Ub the same as Ua.

density of liquid 900 kg/m3
density of solid particle 1100 kg/m3

â–*case 1
kineticTheory off;
g0 =1,

preAlphaExp=500;
packingLimiter on;
expMax 1.0e+03;
alphaMax 0.62;

â–*case 2
kineticTheory on;
g0 = 0;

â–*case 3
kineticTheory off;
g0=0;

As Alberto said, the value of alpha can not be greater than the maximum packing limit,
but the kinetcTheory or particle-particle interaction model did not prevent the void fraction
becoming greate than the maximum packing limit.

What do you suggest?

Thank you in advance.

[hemph]

Hi Guanghao
What is your time step? Packing simulations usually require a very short time step, on the order of 1e-5. For your case 1, note that in OpenFOAM, the particle-particle force has the unit Pa/(kg/m^3), so your g0 should probably be divided by 1100 which might help things. Try to increase the expMax parameter to see an increase in force for an increase in alpha.

//Rasmus