
[Sponsors] 
October 16, 2012, 05:10 
Alpha1 is reducing.....

#1 
Member
Stefano
Join Date: Jul 2009
Posts: 34
Rep Power: 10 
Hi guys,
i have a problem with alpha1 in openFoam 1.6 ext. The value reduces with the time. Here is the output: /**\  =========    \\ / F ield  OpenFOAM Extend Project: Open source CFD   \\ / O peration  Version: 1.6ext   \\ / A nd  Web: www.extendproject.de   \\/ M anipulation   \**/ Build : 1.6ext Exec : interGravitationalFoam Date : Oct 16 2012 Time : 10:56:34 Host : Alderaan PID : 3367 Case : /media/F2D48837D487FFD9/LAVORO/OPENFOAM/REXUS/test nProcs : 1 SigFpe : Enabling floating point exception trapping (FOAM_SIGFPE). // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * // Create time Create mesh for time = 0 Reading g Loading Accelerations.data Loading Omega.data Calculating field g.h Reading inertial Origin Reading field pd Reading field alpha1 Reading field U Reading/calculating face flux field phi Reading transportProperties Selecting incompressible transport model Newtonian Selecting incompressible transport model Newtonian Selecting turbulence model type RASModel Selecting RAS turbulence model kOmegaSST kOmegaSSTCoeffs { alphaK1 0.85034; alphaK2 1; alphaOmega1 0.5; alphaOmega2 0.85616; gamma1 0.5532; gamma2 0.4403; beta1 0.075; beta2 0.0828; betaStar 0.09; a1 0.31; c1 10; } time step continuity errors : sum local = 0, global = 0, cumulative = 0 DICPCG: Solving for pcorr, Initial residual = 0, Final residual = 0, No Iterations 0 time step continuity errors : sum local = 0, global = 0, cumulative = 0 Courant Number mean: 0 max: 0 velocity magnitude: 0 Omega1 initialized to (0 0 0) Omega2 initialized to (0 0 0) Starting time loop Courant Number mean: 0 max: 0 velocity magnitude: 0 deltaT = 0.000595238 Time = 0.000595238 Calculating field g.h Setting accelerations to (112.644 1.09078 0.109779) Setting omega to (0.00875158 0.0815817 0.0820055) Time at instant 1 is: 0 Time at instant 2 is: 0.000595238 Omega at instant 1 is: (0 0 0) Omega at instant 2 is: (0.00875158 0.0815817 0.0820055) Instantaneous angular acceleration is: (14.7027 137.057 137.769) Previous Liquid phase volume fraction = 2.25683e07 Min(alpha1) = 0 Max(alpha1) = 1 DILUPBiCG: Solving for alpha1: solution singularity Liquid phase volume fraction = 2.25683e07 Min(alpha1) = 0 Max(alpha1) = 1 DICPCG: Solving for pd, Initial residual = 1, Final residual = 9.62755e08, No Iterations 329 time step continuity errors : sum local = 1.16224e11, global = 1.58311e23, cumulative = 1.58311e23 time step continuity errors : sum local = 1.16224e11, global = 1.58311e23, cumulative = 3.16623e23 DILUPBiCG: Solving for omega, Initial residual = 0.000655187, Final residual = 1.02048e07, No Iterations 1 DILUPBiCG: Solving for k, Initial residual = 1, Final residual = 1.13822e06, No Iterations 1 ExecutionTime = 27.63 s ClockTime = 27 s Centre of Mass: (3.14552e07 4.44441e05 0.000106692) Courant Number mean: 0.000712776 max: 0.0760438 velocity magnitude: 10.4599 deltaT = 0.000705782 Time = 0.00130102 Calculating field g.h Setting accelerations to (113.489 1.10463 0.123626) Setting omega to (0.00919504 0.0653216 0.0650064) Time at instant 1 is: 0.000595238 Time at instant 2 is: 0.00130102 Omega at instant 1 is: (0.00875158 0.0815817 0.0820055) Omega at instant 2 is: (0.00919504 0.0653216 0.0650064) Instantaneous angular acceleration is: (0.628319 23.0383 24.0855) Previous Liquid phase volume fraction = 2.25683e07 Min(alpha1) = 0 Max(alpha1) = 1 DILUPBiCG: Solving for alpha1, Initial residual = 0.00503263, Final residual = 5.05418e09, No Iterations 2 Liquid phase volume fraction = 2.25683e07 Min(alpha1) = 0 Max(alpha1) = 0.892477 DICPCG: Solving for pd, Initial residual = 0.918833, Final residual = 9.18764e08, No Iterations 329 time step continuity errors : sum local = 1.97884e11, global = 3.14334e22, cumulative = 2.82672e22 time step continuity errors : sum local = 1.97884e11, global = 3.14334e22, cumulative = 5.97006e22 DILUPBiCG: Solving for omega, Initial residual = 0.000743617, Final residual = 2.82979e08, No Iterations 1 DILUPBiCG: Solving for k, Initial residual = 0.333835, Final residual = 2.95525e07, No Iterations 1 ExecutionTime = 46.73 s ClockTime = 47 s Centre of Mass: (1.9406e07 3.56954e05 8.55697e05) Courant Number mean: 0.000676682 max: 0.0735919 velocity magnitude: 8.53713 deltaT = 0.000839638 Time = 0.00214066 Calculating field g.h Setting accelerations to (114.494 1.1211 0.1401) Setting omega to (0.0097226 0.0459778 0.0447833) Time at instant 1 is: 0.00130102 Time at instant 2 is: 0.00214066 Omega at instant 1 is: (0.00919504 0.0653216 0.0650064) Omega at instant 2 is: (0.0097226 0.0459778 0.0447833) Instantaneous angular acceleration is: (0.628319 23.0383 24.0855) Previous Liquid phase volume fraction = 2.25683e07 Min(alpha1) = 0 Max(alpha1) = 0.892477 DILUPBiCG: Solving for alpha1, Initial residual = 0.00309017, Final residual = 2.7438e06, No Iterations 1 Liquid phase volume fraction = 2.25688e07 Min(alpha1) = 1.11844e97 Max(alpha1) = 0.824762 DICPCG: Solving for pd, Initial residual = 0.617329, Final residual = 9.87751e08, No Iterations 287 time step continuity errors : sum local = 1.2306e11, global = 1.7362e24, cumulative = 5.95269e22 time step continuity errors : sum local = 1.2306e11, global = 1.7362e24, cumulative = 5.93533e22 DILUPBiCG: Solving for omega, Initial residual = 0.000872859, Final residual = 9.36405e09, No Iterations 1 DILUPBiCG: Solving for k, Initial residual = 0.263603, Final residual = 1.13085e07, No Iterations 1 ExecutionTime = 64.18 s ClockTime = 64 s Centre of Mass: (1.11854e07 3.07893e05 7.3643e05) Courant Number mean: 0.00056678 max: 0.0512143 velocity magnitude: 5.03397 deltaT = 0.00099707 Time = 0.00313773 Calculating field g.h Setting accelerations to (115.687 1.14066 0.159662) Setting omega to (0.0103491 0.0230069 0.0207683) Time at instant 1 is: 0.00214066 Time at instant 2 is: 0.00313773 Omega at instant 1 is: (0.0097226 0.0459778 0.0447833) Omega at instant 2 is: (0.0103491 0.0230069 0.0207683) Instantaneous angular acceleration is: (0.628319 23.0383 24.0855) Previous Liquid phase volume fraction = 2.25688e07 Min(alpha1) = 1.11844e97 Max(alpha1) = 0.824762 DILUPBiCG: Solving for alpha1, Initial residual = 0.00133358, Final residual = 3.93363e07, No Iterations 1 Liquid phase volume fraction = 2.25688e07 Min(alpha1) = 5.272e08 Max(alpha1) = 0.79509 I have set buoyantPressure for p and pd and I run a modified interFoam. Why do I have this result??? Please help me!!! 

October 16, 2012, 05:16 

#2 
Senior Member
Anton Kidess
Join Date: May 2009
Location: Germany
Posts: 1,261
Rep Power: 23 
You don't give a lot of info, but try using explicit MULES for alpha1. Use a bounded divergence scheme and first order Euler.
__________________
*On twitter @akidTwit *Spend as much time formulating your questions as you expect people to spend on their answer. 

October 16, 2012, 05:36 

#3 
Member
Stefano
Join Date: Jul 2009
Posts: 34
Rep Power: 10 
Which info do you need in particular?
I set MULES as Explicit, and Euler as time scheme. I don't know exactly what you mean with "bounded scheme" Anyway it doesn't change, but still gives back the same error. 

October 16, 2012, 05:42 

#4 
Senior Member
Anton Kidess
Join Date: May 2009
Location: Germany
Posts: 1,261
Rep Power: 23 
Helpful information would be the content of fvSchemes and fvSolution, specifying if you have any special boundary conditions (contact angles, inlets/outlets), or a problematic mesh.
About the bounded schemes  see "4.4.1.1 Schemes for strictly bounded scalar fields" in the user guide. Examples of such schemes are vanLeer01 and Gamma01.  Anton
__________________
*On twitter @akidTwit *Spend as much time formulating your questions as you expect people to spend on their answer. 

October 16, 2012, 05:51 

#5 
Member
Stefano
Join Date: Jul 2009
Posts: 34
Rep Power: 10 
Ok, I attach all the files for boundary condition plus fvSolution and fvScheme.
There is the contact angle acting in the domain. 

October 16, 2012, 06:50 

#6  
Super Moderator
Tobias Holzmann
Join Date: Oct 2010
Location: Leoben (Austria)
Posts: 1,874
Blog Entries: 6
Rep Power: 32 
Quote:
Your Problem description                                           I can not see any problem in your case. You said that the value of alpha1 goes below 1. Thats does not mean to be an error or an numerical problem. I had a look in your time folder 0 and saw that there is no inlet for alpha1 that is set to fixedValue. I have got a notation that ...                                           you are init your alpha1 with the setFields command in a certain region (like the damnBreak tutorial). While solving your case the fraction 1 is distributed in the region (maybe due to gravitation) and your fraction value reaches values below 1 due to distribution and the fact that your fraction 1 is a thin layer or something like that. Further more: You set your BC for U to fixedValue (0 0 0). So I expect that you are in a closed box or sth. like that. Experience                                           A year ago I made a simulation of a chocolate donot. Therefor I set a chocolate field like a plate into my mesh and let it fall down on my donut due to gravitation. I realized that the fraction is distributed due to the gravitation and velocity becouse my chocolate plate was very thin and I used the visco of water and not of chocolate, The results was, that my maxAlpha1 value goes down below 1. That was logical couse the cells which belong to alpha1 = 1 take a high of 4 cells. After the first few timesteps I got values below 1. I had a look at the first calculation steps and realized that problem. To do                                           First of all, write out the first time steps and have a look at it with paraview. It should be possible to see your problem, or if your calculation is correct or wrong. If everything is okay check your boundary conditions, initialization and so on. After that you can go on and check your schemes. Like akidess told you use first order euler schemes for time derivation and bounded scheme for alpha1. But if the problem occure with MULES and DILUPBiCG it seems like the problem I mentioned and that its not a problem in the numerical schemes. I have no idea about your geometry and initialization but I think its due to distribution of the fraction and the value goes under 1 couse your fraction 1 is a thin layer or sth. like that. Further more                                           its possible that I skipped some details?.?.? Without a concrete case I can not help you. I grope on the dark. But I Hope that this will help you a bit Regard Tobi Last edited by Tobi; October 16, 2012 at 07:05. Reason: Adding some lines 

October 16, 2012, 07:04 

#7 
Senior Member
Anton Kidess
Join Date: May 2009
Location: Germany
Posts: 1,261
Rep Power: 23 
I agree with Tobi. I understood your original problem as that the volume integral of alpha1 is reducing (i.e. you are losing mass), but now I realized you meant that your maximum value is reducing (i.e. your field is dispersing). Follow the steps suggested by Tobi.
__________________
*On twitter @akidTwit *Spend as much time formulating your questions as you expect people to spend on their answer. 

October 16, 2012, 10:08 

#9 
Member
Stefano
Join Date: Jul 2009
Posts: 34
Rep Power: 10 
Hi guys,
thanks for your help.... Anyway, the domain is a cylindrical container, with inside some liquid and gas, subjected to acceleration varying with time. Because I have modified a bit the interFoam solver, i send to you my UEqn and pEqn files: maybe i have made some mistakes in writing the code. pEqn { volScalarField rUA = 1.0/UEqn.A(); surfaceScalarField rUAf = fvc::interpolate(rUA); U = rUA*UEqn.H(); surfaceScalarField phiU ( "phiU", (fvc::interpolate(U) & mesh.Sf()) + fvc::ddtPhiCorr(rUA, rho, U, phi) ); adjustPhi(phiU, U, pd); phi = phiU + ( fvc::interpolate(interface.sigmaK())*fvc::snGrad(a lpha1)  ghf*fvc::snGrad(rho) )*rUAf*mesh.magSf(); for(int nonOrth=0; nonOrth<=nNonOrthCorr; nonOrth++) { fvScalarMatrix pdEqn ( fvm::laplacian(rUAf, pd) == fvc::div(phi) ); pdEqn.setReference(pdRefCell, pdRefValue); if (corr == nCorr1 && nonOrth == nNonOrthCorr) { pdEqn.solve(mesh.solver(pd.name() + "Final")); } else { pdEqn.solve(mesh.solver(pd.name())); } if (nonOrth == nNonOrthCorr) { phi = pdEqn.flux(); } } U += rUA*fvc::reconstruct((phi  phiU)/rUAf); U.correctBoundaryConditions(); # include "continuityErrs.H" p = pd + rho*gh; if (pd.needReference()) { p += dimensionedScalar ( "p", p.dimensions(), pRefValue  getRefCellValue(p, pdRefCell) ); } } UEqn surfaceScalarField muEff ( "muEff", twoPhaseProperties.muf() + fvc::interpolate(rho*turbulence>nut()) ); fvVectorMatrix UEqn ( fvm::ddt(rho, U) + fvm::div(rhoPhi, U)  fvm::laplacian(muEff, U)  (fvc::grad(U) & fvc::grad(muEff)) + Fcentrifugal + Feuler + Fcoriolis // fvc::div(muEff*(fvc::interpolate(dev(fvc::grad(U)) ) & mesh.Sf())) ); UEqn.relax(); if (momentumPredictor) { solve ( UEqn == fvc::reconstruct ( ( fvc::interpolate(interface.sigmaK())*fvc::snGrad(a lpha1)  ghf*fvc::snGrad(rho)  fvc::snGrad(pd) ) * mesh.magSf() ) ); } I have also written the following lines in the main code, to take the gravity variation in account: dimensionedVector acc("Acceleration", dimAcceleration,it(runTime.timeOutputValue())); g = acc; Info<< "Calculating field g.h\n" << endl; volScalarField gh("gh", g & mesh.C()); surfaceScalarField ghf("gh", g & mesh.Cf()); Maybe now my problem is clearer (I hope). 

Tags 
alpha1 
Thread Tools  
Display Modes  


Similar Threads  
Thread  Thread Starter  Forum  Replies  Last Post 
Formulation in compressibleInterFoam  scttmllr  OpenFOAM Running, Solving & CFD  57  August 25, 2017 06:03 
Multiple floating objects  CKH  OpenFOAM Running, Solving & CFD  12  March 21, 2016 14:05 
interFoam VOF is loosing fluid  wersoe  OpenFOAM Running, Solving & CFD  12  June 26, 2013 08:13 
Problem with interFoam; Wave/wiggle alpha1 behavior  JonW  OpenFOAM  3  February 23, 2013 21:41 
Fixing the alpha1 min/max under/overshoots in interFoam.C solver.  idrama  OpenFOAM  1  January 27, 2010 20:34 