CFD Online Logo CFD Online URL
www.cfd-online.com
[Sponsors]
Home > Forums > OpenFOAM Programming & Development

Manipulating the alpha1-field after every time step

Register Blogs Members List Search Today's Posts Mark Forums Read

Reply
 
LinkBack Thread Tools Display Modes
Old   May 8, 2013, 03:36
Default Manipulating the alpha1-field after every time step
  #1
Member
 
Join Date: Aug 2011
Posts: 74
Rep Power: 5
idefix is on a distinguished road
Hello,

I am using the interFoam-solver. Alpha1 is characterizing the volume fraction (alpha1= 0 means: there is only gas in the cell, alpha1=1 means: there is only liquid in the cell).
Now Ive got the problem, that alpha1 is getting negative or greater than 1 which is not realistic.

My idea is to "correct" the alpha1 after every time step.
I need something which goes into every cell: If alpha1 is negative or nearly 0 (for example 0.1e-6), alpha1 should be put to 0. If alpha1 is greater than or nearly 1, alpha1 should be put to 1.

I would put this function in tha alphaEn.H - file after the for loop.

Can anybody help me?

Thanks a lot
idefix is offline   Reply With Quote

Old   May 8, 2013, 04:02
Default
  #2
Senior Member
 
Pablo Higuera
Join Date: Jan 2011
Posts: 237
Rep Power: 7
Phicau is on a distinguished road
Hi,

I would not suggest doing this, as MULES handles itself the boundedness, continuity and conservation of alpha1. Although you may obtain vales slightly lower than 0 or larger than 1, if your BCs are correctly set, these will only be on the order of machine precision.

However, to do what you asked:

Code:
threshold = 0.01; // You name it
alpha1 = alpha1 * pos(alpha1-threshold) * neg(alpha1-1+threshold) + pos(alpha1-1+threshold)
Best,

Pablo
Phicau is offline   Reply With Quote

Old   May 13, 2013, 11:03
Default
  #3
Member
 
Join Date: Nov 2012
Location: Liverpool, UK
Posts: 51
Rep Power: 4
sophie_l is on a distinguished road
Hi Pablo,

I am using interFoam/waveFoam as well. I remember in the code, the value of alpha1 is restricted between 0 and 1. However, the case I am running at the moment gives out results like this, where Max(alpha1) = 1.57337. But the case is still running even though.

Quote:
Courant Number mean: 0.00685244 max: 0.558313
Interface Courant Number mean: 6.62914e-05 max: 0.558313
deltaT = 0.00132421
Time = 0.186161

MULES: Solving for alpha1
Liquid phase volume fraction = 0.666668 Min(alpha1) = 0 Max(alpha1) = 1.16009
MULES: Solving for alpha1
Liquid phase volume fraction = 0.666668 Min(alpha1) = 0 Max(alpha1) = 1.57337
GAMG: Solving for p_rgh, Initial residual = 5.1474e-05, Final residual = 6.30767e-06, No Iterations 1
time step continuity errors : sum local = 1.13037e-05, global = -3.39547e-08, cumulative = -6.81206e-06
DILUPBiCG: Solving for epsilon, Initial residual = 0.00528157, Final residual = 7.39422e-09, No Iterations 7
DILUPBiCG: Solving for k, Initial residual = 0.00656354, Final residual = 3.16796e-09, No Iterations 8
ExecutionTime = 1408.43 s ClockTime = 1465 s
I am thinking it should be due to the input boundary conditions, however, I haven't got a clue where the check should start with. Do you some suggestions on this please?

Thanks a lot!
Sophie
sophie_l is offline   Reply With Quote

Old   May 13, 2013, 11:30
Default
  #4
Senior Member
 
Pablo Higuera
Join Date: Jan 2011
Posts: 237
Rep Power: 7
Phicau is on a distinguished road
Hi

it is very likely that your alpha1 unboundedness is caused by the BCs. I suggest you take a look at this post for BC reference:

Waves2Foam Related Topics

if that does not solve your problems, explain your case there and I am sure they will be able to help you.

Best,

Pablo
Phicau is offline   Reply With Quote

Old   May 13, 2013, 15:24
Default
  #5
ngj
Senior Member
 
Niels Gjoel Jacobsen
Join Date: Mar 2009
Location: Deltares, Delft, The Netherlands
Posts: 1,607
Rep Power: 25
ngj will become famous soon enoughngj will become famous soon enough
Hi Sophie,

From the snap-shot of the log file, it is clear that:

1. You are not solving the momentum equation.

2. You are only solving the pressure equation one (1) time per time step.

If you relate that with the fact that interFoam/waveFoam is based on a PIMPLE-algorithm, then you will have no chance what so ever in obtaining mass conservation on the cell level; in your case the error in mass conservation is considerable. This is most probably the reason for the unboundedness of alpha1.

So: start solving for the momentum equation and have 3-4 iterations in the pressure equation per time step. This will greatly improve on your mass conservation problems and therefore improve on the boundedness of alpha1.

Kind regards,

Niels
ngj is offline   Reply With Quote

Old   May 13, 2013, 16:00
Default
  #6
Member
 
Join Date: Nov 2012
Location: Liverpool, UK
Posts: 51
Rep Power: 4
sophie_l is on a distinguished road
Hi Niels,

Thank you so much for your reply. Really appreciate your help so far!

That's quite inspiring really. I'll resort to the code to work it out. By the way, I am using OF 1.7.1 based on PISO algorithm. However, in the cases I ran earlier, I didn't see alpha1 could be larger than 1.5, at most 1.0*.

Many thanks,
Sophie
sophie_l is offline   Reply With Quote

Old   May 13, 2013, 16:14
Default
  #7
ngj
Senior Member
 
Niels Gjoel Jacobsen
Join Date: Mar 2009
Location: Deltares, Delft, The Netherlands
Posts: 1,607
Rep Power: 25
ngj will become famous soon enoughngj will become famous soon enough
Hi Sophie,

If you are using the PISO algorithm, the problem is at least as bad as using PIMPLE (though, if your relaxing factors are 1, then the results should be identical).

You do not need to modify/read the code, but merely change some keywords in fvSolution.

Furthermore, I can recommend that you read the section on the PISO algortihm in

Code:
@book { Ferziger2002,
Author = {Ferziger, J. H. and Peric, M},
Title = {Computational Methods for Fluid Dynamics},
Publisher = {Springer},
Year = {2002},
Edition = {3rd},
}
which will aid you in understand the importance of the settings for the PISO-algorithm.

Kind regards,

Niels
ngj is offline   Reply With Quote

Old   May 13, 2013, 16:27
Default
  #8
Member
 
Join Date: Nov 2012
Location: Liverpool, UK
Posts: 51
Rep Power: 4
sophie_l is on a distinguished road
Hi Niels,

Thanks a lot!

Fortunately, I've borrowed this book from the library. I'll read this part right now.

Cheers,
Sophie
sophie_l is offline   Reply With Quote

Old   May 14, 2013, 03:01
Default
  #9
Member
 
Join Date: Aug 2011
Posts: 74
Rep Power: 5
idefix is on a distinguished road
Hi,

the discussion is really helpful for me and I have got one question to it.

My simulations are not running very well. Soon after the beginning the velocities are increasing and as a result my time step is decreasing till the simulation stopps ( I am using adjustableRunTime). The turbulent viscosity is getting very small too.
I am using the k-epsilon-model.
At the beginning I thought that alpha1 is the reason for this behaviour but the discussion here gives some hints that the problem could be somewhere else.

Is there a possiblity to change something in the fvSolution-file so that the simulation get stable?

At the moment my fvSolution-file looks like:
Quote:
solvers
{
pcorr
{
solver PCG;
preconditioner DIC;
tolerance 1e-10;
relTol 0;
}

p_rgh
{
solver PCG;
preconditioner DIC;
tolerance 1e-07;
relTol 0.05;
}

p_rghFinal
{
solver PCG;
preconditioner DIC;
tolerance 1e-07;
relTol 0;
}

"(U|k|epsilon)"
{
solver PBiCG;
preconditioner DILU;
tolerance 1e-06;
relTol 0;
}

"(U|k|epsilon)Final"
{
solver PBiCG;
preconditioner DILU;
tolerance 1e-08;
relTol 0;
}
}
PIMPLE
{
momentumPredictor no;
nCorrectors 3;
nNonOrthogonalCorrectors 0;
nAlphaCorr 1;
nAlphaSubCycles 4;
cAlpha 2;
}
My log-file is:
Quote:
Courant Number mean: 0.001746147011 max: 0.512101182015
Interface Courant Number mean: 2.50102364228e-05 max: 0.0844578543383
deltaT = 3.4147974996e-18
Time = 0.0041231153508954621

MULES: Solving for alpha1
Phase-1 volume fraction = 0.0910985755396 Min(alpha1) = -8.28923892037e-44 Max(alpha1) = 1.00000001104
MULES: Solving for alpha1
Phase-1 volume fraction = 0.0910985743502 Min(alpha1) = -3.51986718636e-44 Max(alpha1) = 1.00000001103
MULES: Solving for alpha1
Phase-1 volume fraction = 0.091098573153 Min(alpha1) = -4.940549133e-45 Max(alpha1) = 1.00000001103
MULES: Solving for alpha1
Phase-1 volume fraction = 0.0910985719483 Min(alpha1) = -3.27296438965e-44 Max(alpha1) = 1.00000001102
DICPCG: Solving for p_rgh, Initial residual = 0.0269677660092, Final residual = 0.00113281797768, No Iterations 2
time step continuity errors : sum local = 3.66942289335e-07, global = 1.44616781464e-07, cumulative = 0.000266660876022
DICPCG: Solving for p_rgh, Initial residual = 0.00387752776814, Final residual = 0.000186226585041, No Iterations 13
time step continuity errors : sum local = 6.10538629187e-08, global = 1.86907321117e-08, cumulative = 0.000266679566754
DICPCG: Solving for p_rgh, Initial residual = 0.000664946844126, Final residual = 9.53727056345e-08, No Iterations 120
time step continuity errors : sum local = 3.12780800404e-11, global = 1.96037099123e-12, cumulative = 0.000266679568714
DILUPBiCG: Solving for epsilon, Initial residual = 0.00538752754077, Final residual = 4.08274384377e-10, No Iterations 2
DILUPBiCG: Solving for k, Initial residual = 0.00378438227684, Final residual = 2.31895835788e-09, No Iterations 1
max(k) = max(k) [0 2 -2 0 0 0 0] 2.60958163608e+13 min(k) = min(k) [0 2 -2 0 0 0 0] 0.00172278530058
max(epsilon) = max(epsilon) [0 2 -3 0 0 0 0] 2.84510295356e+22 min(epsilon) = min(epsilon) [0 2 -3 0 0 0 0] 15.772024055
max(U) = max(U) [0 1 -1 0 0 0 0] (1.3707489336e+12 2.3034155804e+12 9.11378593063e+12) min(U) = min(U) [0 1 -1 0 0 0 0] (-1.87325084502e+12 -3.151035603e+12 -5.6941716128e+12

ExecutionTime = 29474.2 s ClockTime = 29593 s
Thanks a lot
idefix is offline   Reply With Quote

Old   May 14, 2013, 03:55
Default
  #10
ngj
Senior Member
 
Niels Gjoel Jacobsen
Join Date: Mar 2009
Location: Deltares, Delft, The Netherlands
Posts: 1,607
Rep Power: 25
ngj will become famous soon enoughngj will become famous soon enough
Hi Idefix,

Try running your simulation without a turbulence closure. Your turbulence quantities are diverged, so poor boundary conditions on k/epsilon could be the source of your crash.

Kind regards

Niels
ngj is offline   Reply With Quote

Old   May 14, 2013, 04:29
Default
  #11
Member
 
Join Date: Aug 2011
Posts: 74
Rep Power: 5
idefix is on a distinguished road
Hi Niels,

thanks for your answer but I am not so familiar with OpenFoam.
How can I run my simulation without turbulence closure?
Going away from URANS?

thanks again
Idefix
idefix is offline   Reply With Quote

Old   May 14, 2013, 05:07
Default
  #12
ngj
Senior Member
 
Niels Gjoel Jacobsen
Join Date: Mar 2009
Location: Deltares, Delft, The Netherlands
Posts: 1,607
Rep Power: 25
ngj will become famous soon enoughngj will become famous soon enough
Hi Idefix,

Open the file

Code:
constant/RASProperties
and change your current RASModel to

Code:
RASModel    laminar;
Kind regards

Niels
ngj is offline   Reply With Quote

Old   May 14, 2013, 06:34
Default
  #13
Member
 
Join Date: Aug 2011
Posts: 74
Rep Power: 5
idefix is on a distinguished road
Hi Niels,

in my case Ive got a pipe flow and the Reynolds-numbers says that I have got a turbulent flow.

For my understanding I should calculate the flow turbulent.

Is there something else I can try?

Thanks a lot,
Idefix
idefix is offline   Reply With Quote

Old   May 14, 2013, 06:53
Default
  #14
ngj
Senior Member
 
Niels Gjoel Jacobsen
Join Date: Mar 2009
Location: Deltares, Delft, The Netherlands
Posts: 1,607
Rep Power: 25
ngj will become famous soon enoughngj will become famous soon enough
Hi Idefix,

My suggestion is only meant for debugging. Turning off the turbulence, you can lower the complexity and see, whether the turbulence model is the cause of your problems. If yes, you know that you need to change e.g. the boundary conditions to stabilise the model.

You could also try to use upwind on all of the divergence schemes in fvSchemes.

Kind regards

Niels
ngj is offline   Reply With Quote

Reply

Thread Tools
Display Modes

Posting Rules
You may not post new threads
You may not post replies
You may not post attachments
You may not edit your posts

BB code is On
Smilies are On
[IMG] code is On
HTML code is Off
Trackbacks are On
Pingbacks are On
Refbacks are On


Similar Threads
Thread Thread Starter Forum Replies Last Post
mixerVesselAMI2D's mass is not balancing sharonyue OpenFOAM Running, Solving & CFD 6 June 10, 2013 09:34
Orifice Plate with a fully developed flow - Problems with convergence jonmec OpenFOAM Running, Solving & CFD 3 July 28, 2011 05:24
Extrusion with OpenFoam problem No. Iterations 0 Lord Kelvin OpenFOAM 6 April 12, 2011 11:24
Full pipe 3D using icoFoam cyberbrain OpenFOAM 4 March 16, 2011 10:20
Differences between serial and parallel runs carsten OpenFOAM Bugs 11 September 12, 2008 11:16


All times are GMT -4. The time now is 05:55.