# alpha in twoPhaseEulerFoam

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

 March 15, 2016, 19:44 alpha in twoPhaseEulerFoam #1 Senior Member   Join Date: Jan 2013 Posts: 368 Rep Power: 10 Hi, IN twoPhaseEulerFoam of OF300 version, we have following lines in createFields.H: Code: ```phaseModel& phase1 = fluid.phase1(); phaseModel& phase2 = fluid.phase2(); volScalarField& alpha1 = phase1; volScalarField& alpha2 = phase2;``` I think alpha1 and alpha2 correspond to the volume fraction of dispersed and continuous phases. Here phase1 and phase2 are the object of class phaseModel, then why can we still specify them to alpha1 and alpha2, which are the objects of class volScalarField? Confused here. Thank you.

 March 15, 2016, 20:16 #2 Senior Member   Join Date: Jan 2013 Posts: 368 Rep Power: 10 Hi, Another question about the twoPhaseEulerFoam is for the quantity residualAlpha: Code: ``` //- Return the residual phase-fraction for given phase // Used to stabilize the phase momentum as the phase-fraction -> 0 const dimensionedScalar& residualAlpha() const { return residualAlpha_; }``` So in the computations of particle laden isothermal jet flows, how low is fine typically? The maximum value of alpha_particle is always at the jet exit, say 1.0e-5. Thank you.

March 16, 2016, 04:23
#3
Member

Juho Peltola
Join Date: Mar 2009
Location: Finland
Posts: 89
Rep Power: 13
Quote:
 Originally Posted by openfoammaofnepo Hi, IN twoPhaseEulerFoam of OF300 version, we have following lines in createFields.H: Code: ```phaseModel& phase1 = fluid.phase1(); phaseModel& phase2 = fluid.phase2(); volScalarField& alpha1 = phase1; volScalarField& alpha2 = phase2;``` I think alpha1 and alpha2 correspond to the volume fraction of dispersed and continuous phases. Here phase1 and phase2 are the object of class phaseModel, then why can we still specify them to alpha1 and alpha2, which are the objects of class volScalarField? Confused here. Thank you.
phaseModel is derived from volScalarField and it is the volume fraction field.
https://github.com/OpenFOAM/OpenFOAM...l/phaseModel.H

Quote:
 Originally Posted by openfoammaofnepo Another question about the twoPhaseEulerFoam is for the quantity residualAlpha. So in the computations of particle laden isothermal jet flows, how low is fine typically? The maximum value of alpha_particle is always at the jet exit, say 1.0e-5. Thank you.
It depends and you will have to test it yourself. You can keep lowering the residual value until it is low enough for you or numerical anomalies start to appear. At very low volume fractions Lagrangian tracking should also be an option?

Probably not an issue for an isothermal jet, but if compressibility effects are significant continuity errors may also appear when volume fraction drops below the dgdt residual alphas that are hard coded (1e-4) into the twoPhaseSystem:

https://github.com/OpenFOAM/OpenFOAM...oPhaseSystem.C

Code:
```        forAll(dgdt_, celli)
{
if (dgdt_[celli] > 0.0)
{
Sp[celli] -= dgdt_[celli]/max(1.0 - alpha1[celli], 1e-4);
Su[celli] += dgdt_[celli]/max(1.0 - alpha1[celli], 1e-4);
}
else if (dgdt_[celli] < 0.0)
{
Sp[celli] += dgdt_[celli]/max(alpha1[celli], 1e-4);
}
}```

 March 16, 2016, 04:24 #4 Senior Member   Gerhard Holzinger Join Date: Feb 2012 Location: Austria Posts: 287 Rep Power: 23 If you have a look at the file phaseModel.H (e.g. here: https://github.com/OpenFOAM/OpenFOAM...l/phaseModel.H), you will see that the class phaseModel is derived (among others) from volScalarField. Thus, the class phaseModel is the volume fraction field alpha. See also the constructor in phaseModel.C for a further clue.

 March 16, 2016, 05:55 #5 Senior Member   Join Date: Jan 2013 Posts: 368 Rep Power: 10 Dear Gerhard and Juho, Thank you for your help. Make sense! For the turbulence modelling in twoPhaseEulerFoam, we have the following: Code: ` if (pimple.turbCorr()) { fluid.correctTurbulence(); }` In TwoPhaseSystem.C, we have see this: Code: `void Foam::twoPhaseSystem::correctTurbulence() { phase1_.turbulence().correct(); phase2_.turbulence().correct(); }` Phase1_. and phase2_ are the objects of class phaseModel. So in phaseModel.C, we have : Code: `Foam::PhaseCompressibleTurbulenceModel& Foam::phaseModel::turbulence() { return turbulence_(); }` Since turbulence_ is the pointer to the PhaseCompressibleTurbulenceModel. However, in the class of PhaseCompressibleTurbulenceModel, I cannot see further what is the function of correct(). Could you please help to point out where it is defined? Thank you so much.

 March 16, 2016, 06:20 #6 Senior Member   Gerhard Holzinger Join Date: Feb 2012 Location: Austria Posts: 287 Rep Power: 23 phase2_.turbulence().correct(); the method turbulence() returns the pointer to the turbulence model. And on this pointer the method correct() is called. The turbulence models used for the phases are of the type phaseCompressibleTurbulenceModel, and here the interesting part begins. It is true that you will not find any method correct() in phaseCompressibleTurbulenceModel. This is because this turbulence model class is derived from other base classes. The templated turbulence model framework is designed to cover everything from single to multiphase and from incompressible to compressible flows. In order to find the method correct() you have been looking so far, you need to dig through the organisation of the turbulence model framework. Maybe Section 20 on turbulence models is of help https://github.com/ParticulateFlow/O...Manual_PFM.pdf

March 16, 2016, 07:10
#7
Senior Member

Join Date: Jan 2013
Posts: 368
Rep Power: 10
Thank you Gerhard!

cheer,
OFFO

Quote:
 Originally Posted by GerhardHolzinger phase2_.turbulence().correct(); the method turbulence() returns the pointer to the turbulence model. And on this pointer the method correct() is called. The turbulence models used for the phases are of the type phaseCompressibleTurbulenceModel, and here the interesting part begins. It is true that you will not find any method correct() in phaseCompressibleTurbulenceModel. This is because this turbulence model class is derived from other base classes. The templated turbulence model framework is designed to cover everything from single to multiphase and from incompressible to compressible flows. In order to find the method correct() you have been looking so far, you need to dig through the organisation of the turbulence model framework. Maybe Section 20 on turbulence models is of help https://github.com/ParticulateFlow/O...Manual_PFM.pdf

March 16, 2016, 18:21
#8
Senior Member

Join Date: Jan 2013
Posts: 368
Rep Power: 10
Dear Gerhard,

I spent some time in reading your document. This is the best OF help file I ever read, except the OF theory guide and user's guide.

Still about the turbulence model, I can understand what you mean about the implementations of turbulence model in twoPhaseEulerFoam. However, I think it is a little difficult to see the role of the souece file: phaseCompressibleTurbulenceModels.C.

I think it is a connection between the current solver and the turbulence model library in /src. Could you please say more about how they are linked? Thank you.

Quote:
 Originally Posted by GerhardHolzinger phase2_.turbulence().correct(); the method turbulence() returns the pointer to the turbulence model. And on this pointer the method correct() is called. The turbulence models used for the phases are of the type phaseCompressibleTurbulenceModel, and here the interesting part begins. It is true that you will not find any method correct() in phaseCompressibleTurbulenceModel. This is because this turbulence model class is derived from other base classes. The templated turbulence model framework is designed to cover everything from single to multiphase and from incompressible to compressible flows. In order to find the method correct() you have been looking so far, you need to dig through the organisation of the turbulence model framework. Maybe Section 20 on turbulence models is of help https://github.com/ParticulateFlow/O...Manual_PFM.pdf

March 16, 2016, 19:01
#9
Senior Member

Join Date: Jan 2013
Posts: 368
Rep Power: 10
Hello,

About the thermophysical model, let us say particle laden flows, in the tutorials of fluidizedbed in twophaseEulerFoam (in OF3.0.0). the Theromotype for air (continuous phase) is
Code:
```thermoType
{
type            heRhoThermo;
mixture         pureMixture;
transport       const;
thermo          hConst;
equationOfState perfectGas;
specie          specie;
energy          sensibleInternalEnergy;
}```
For particle, it is
Code:
```thermoType
{
type            heRhoThermo;
mixture         pureMixture;
transport       const;
thermo          hConst;
equationOfState rhoConst;
specie          specie;
energy          sensibleInternalEnergy;
}```
So we can see that for particle density (here I think it is talking about the material density for a single particle), it is unchanged during the simulation. However, for the twoPhaseEulerFoam (OF3.0.0), the energy equation is solved and we have the following to update the thermo-related parameters:

Code:
```    thermo1.correct();
Info<< "min " << thermo1.T().name()
<< " " << min(thermo1.T()).value() << endl;

thermo2.correct();
Info<< "min " << thermo2.T().name()
<< " " << min(thermo2.T()).value() << endl;```
Since the thermoType is heRhoThermo, so in heRhoThermo.C, we can see what quantities are updated:
Code:
```  scalarField& TCells = this->T_.internalField();
scalarField& psiCells = this->psi_.internalField();
scalarField& rhoCells = this->rho_.internalField();
3>    scalarField& muCells = this->mu_.internalField();
scalarField& alphaCells = this->alpha_.internalField();```
My questions are:
(1) all these quantities are field quantities, how are they linked with the calculations for particulars. For example, rho1 are the dispersed phase meterial density.
(2) I think rho in heRhoThermo.C is updated with EoS with the following:
Code:
`rhoCells[celli] = mixture_.rho(pCells[celli], TCells[celli]);`
However, the the fluidizedbed case file, it is rhoConst. So why can we still use
Code:
`thermo1.correct(); // update the thermophysical parameters of dispersed phase`
to update it.

Confused. Thank you!

Quote:
 Originally Posted by juho phaseModel is derived from volScalarField and it is the volume fraction field. https://github.com/OpenFOAM/OpenFOAM...l/phaseModel.H It depends and you will have to test it yourself. You can keep lowering the residual value until it is low enough for you or numerical anomalies start to appear. At very low volume fractions Lagrangian tracking should also be an option? Probably not an issue for an isothermal jet, but if compressibility effects are significant continuity errors may also appear when volume fraction drops below the dgdt residual alphas that are hard coded (1e-4) into the twoPhaseSystem: https://github.com/OpenFOAM/OpenFOAM...oPhaseSystem.C Code: ``` forAll(dgdt_, celli) { if (dgdt_[celli] > 0.0) { Sp[celli] -= dgdt_[celli]/max(1.0 - alpha1[celli], 1e-4); Su[celli] += dgdt_[celli]/max(1.0 - alpha1[celli], 1e-4); } else if (dgdt_[celli] < 0.0) { Sp[celli] += dgdt_[celli]/max(alpha1[celli], 1e-4); } }```

 March 17, 2016, 04:11 #10 Senior Member   Gerhard Holzinger Join Date: Feb 2012 Location: Austria Posts: 287 Rep Power: 23 Hi, I am not sure whether I understand you question correctly. The framework of OpenFOAM treats certain quantities as if they were variables, e.g. the phase density. The framework then provides methods to evaluate these variables, e.g. the equation of state for phase-densities. Finally, the framework offers a wide choice of models, which range from trivial (e.g. rhoConst) to more comples (e.g. perfectGas). You can see this all over OpenFOAM, e.g. the liftModel framework offers a model named noLift, which returns a field of zeros for the lift coefficient. Thus, the lift force evaluates to zero. It is a little wasted effort to multiply stuff by zero to get zero, however, this keeps the framework simple and versatile. There are no if conditions related to special cases such as laminar flow, no lift force or constant phase density. So, don't be confused when you see a method call to correct(), solve() or some other evaluation, when there is nothing to evaluate. If you go down the sequence of function calls, with the apropriate classes, then you will probably find methods which do nothing, in cases when it is appropriate to do nothing. openfoammaofnepo and jylee4 like this.

 March 18, 2016, 20:03 #11 Senior Member   Join Date: Jan 2013 Posts: 368 Rep Power: 10 Thank you so much for your help. Based on the tutorials, I can see that for particle-laden flows, different turbulence models are used for these two phases. For example, the available LES model include: Code: ```//...List of LES models // Smagorinsky // kEqn // SmagorinskyZhang // NicenoKEqn // continuousGasKEqn // kineticTheoryModel // phasePressureModel``` Do you have any comments or suggestions for the selection of these models for the two phases and their performance? Can I use Smagorinsky for both phases? --------------- The references are listed below for the future checking: 1, NicenoKEqn: http://www.sciencedirect.com/science...0925090800208X No turbulence model is used (assume laminar) for the dispersed phase (gas phase), while for continuous phase (liquid) two models (dynamics smagorinsky model and one equation model) are discussed. 2, SmagorinskyZhang http://www.sciencedirect.com/science...09250906005446 This belongs to a modified version of the Smagorinski model. The bubble induced turbulence is introduced into the effective sub-grid scale viscosity. Meanwhile, the effective viscosity of dispersed phase is estimated from that of the continuous phase with local density ratio. lourencosm likes this. Last edited by openfoammaofnepo; March 20, 2016 at 18:00.

 March 19, 2016, 05:27 #12 Senior Member   Gerhard Holzinger Join Date: Feb 2012 Location: Austria Posts: 287 Rep Power: 23 Note the difference between which models OpenFOAM offers you to use, and what models make sense for a particular application. For the choice of models you best consult relevant academic literature for your particular application. twoPhaseEulerFoam is completely agnostic of the actual physical experiment/set-up you want to simulate. Nothing in OpenFOAM prevents you from applying the kineticTheoryModel to bubbly flows or the SmagorinskyZhang model for gas-solid flows. In the case set-up this selection would be allowed, however, the kineticTheoryModel was derived specifically for gas-solid flows and the SmagorinskyZhang for bubbly flows. openfoammaofnepo, lourencosm, jylee4 and 1 others like this.

March 20, 2016, 18:19
#13
Senior Member

Join Date: Jan 2013
Posts: 368
Rep Power: 10
Based on my investigations, it seems that there is no suitable LES model in OF so far when the dispersed phase is solid particle, such as in particle laden flows. This is based on the latest version of OpenFOAM 3.0.1. The exisiting LES models seems more suitable for gas-liquid flows, such as bubbly flows.

Quote:
 Originally Posted by GerhardHolzinger Note the difference between which models OpenFOAM offers you to use, and what models make sense for a particular application. For the choice of models you best consult relevant academic literature for your particular application. twoPhaseEulerFoam is completely agnostic of the actual physical experiment/set-up you want to simulate. Nothing in OpenFOAM prevents you from applying the kineticTheoryModel to bubbly flows or the SmagorinskyZhang model for gas-solid flows. In the case set-up this selection would be allowed, however, the kineticTheoryModel was derived specifically for gas-solid flows and the SmagorinskyZhang for bubbly flows.

 March 21, 2016, 10:32 #14 Senior Member   Join Date: Jan 2013 Posts: 368 Rep Power: 10 Hi, In twoPhaseEulerFoam, is it possible to print out the residual and number of the iterations for the momentum equations as well? In stead of just say: 'Constructing momentum equations'. Since I should know the residuals of solving the momentum equations. In the current solver of twoPhaseEulerFoam, we can get as follows: ========== PIMPLE: iteration 1 MULES: Solving for alpha.particles MULES: Solving for alpha.particles alpha.particles volume fraction = 0 Min(alpha.particles) = 0 Max(alpha.particles) = 2.4e-05 Constructing momentum equations min T.particles 300 min T.air 300 GAMG: Solving for p_rgh, Initial residual = 1, Final residual = 6.98034e-09, No Iterations 22 ========== Thank you so much.

March 22, 2016, 07:39
#15
Member

Join Date: May 2015
Posts: 34
Rep Power: 7
Quote:
 Originally Posted by openfoammaofnepo Hi, In twoPhaseEulerFoam, is it possible to print out the residual and number of the iterations for the momentum equations as well?
Hi,
I am not 100% sure but, I think the only way this might be possible is with the debug flag in the compilation of openfoam. In pimpleControl.C is the following which indicates printing of more variables with their residuals:
Code:
```if (debug)
{
Info<< algorithmName_ << " loop:" << endl;

Info<< "    " << variableName
<< " PIMPLE iter " << corr_
<< ": ini res = "
<< residualControl_[fieldI].initialResidual
<< ", abs tol = " << residual
<< " (" << residualControl_[fieldI].absTol << ")"
<< ", rel tol = " << relative
<< " (" << residualControl_[fieldI].relTol << ")"
<< endl;
}```
However, I haven't tried this, and I am not sure if this will give velocity residuals, as velocity is not solved but corrected in the solution scheme. Maybe someone has tested this or has a better answer.

March 22, 2016, 08:23
#16
Senior Member

Gerhard Holzinger
Join Date: Feb 2012
Location: Austria
Posts: 287
Rep Power: 23
What you describe has nothing to do with debug flags in compilation.
The code you mentioned is part of OpenFOAM's run-time debugging mechanism.

It is enabled by adding the following dictionary to controlDict. Check out the global controlDict for other available options.
The DebugSwitches dictionary in your case's controlDict overrides the settings of the global controlDict. Thus, there is no need to edit the global controlDict. Unless you know what you are doing, and have ample reason to change OpenFOAM's standard behaviour.

Code:
```DebugSwitches
{
pimpleControl        1;
}```

Quote:
 Originally Posted by mnikku Hi, I am not 100% sure but, I think the only way this might be possible is with the debug flag in the compilation of openfoam. In pimpleControl.C is the following which indicates printing of more variables with their residuals: Code: ```if (debug) { Info<< algorithmName_ << " loop:" << endl; Info<< " " << variableName << " PIMPLE iter " << corr_ << ": ini res = " << residualControl_[fieldI].initialResidual << ", abs tol = " << residual << " (" << residualControl_[fieldI].absTol << ")" << ", rel tol = " << relative << " (" << residualControl_[fieldI].relTol << ")" << endl; }```

Furthermore, twoPhaseEulerFoam does not solve for the momentum explicitely. Thus, you won't find any calls similar to solve(UEqn).
The solution algorithm is described in Henrik Rusche's thesis [http://powerlab.fsb.hr/ped/kturbo/Op...chePhD2002.pdf].

 May 3, 2016, 11:25 #17 New Member   Ehsan Join Date: Aug 2010 Location: QC, Canada Posts: 29 Rep Power: 12 Hi friends, Regarding to existence of very low value of alpha phase fraction (0.88e-3) , I was recently noticed twoPhaseEulerFoam doesnt solve alpha equation . I think I should change this section somehow (twoPhaseSystem.C) to enable code see alpha: Code: ``` forAll(dgdt_, celli) { if (dgdt_[celli] > 0.0 && alpha1[celli] > 0.0) { Sp[celli] -= dgdt_[celli]*alpha1[celli]; Su[celli] += dgdt_[celli]*alpha1[celli]; } else if (dgdt_[celli] < 0.0 && alpha1[celli] < 1.0) { Sp[celli] += dgdt_[celli]*(1.0 - alpha1[celli]); } }``` I would be appreciated if u give me your comments. P.S I am working on an electrochemical reactor (please see attached pics). Cylinder vessel is occupied completely by water. Inlet is wallRotatingvelocity type and gas phase fraction (air) is equal to 0.88e-3 in Inlet. By the way, bubble size diameter is 145e-6 (order of micro) . gemometry.jpg

 May 3, 2016, 12:02 #18 Senior Member   Gerhard Holzinger Join Date: Feb 2012 Location: Austria Posts: 287 Rep Power: 23 The method twoPhaseSystem::solve() contains the solution of the alpha equation. Like in this section: Code: ``` MULES::explicitSolve ( geometricOneField(), alpha1, phi_, alphaPhic1, Sp, Su, phase1_.alphaMax(), 0 );``` The solution of the transport equation for the volume fraction looks now a bit more compilicated than in former days, see https://github.com/OpenFOAM/OpenFOAM...oam/alphaEqn.H However, the volume fraction is still solved for, otherwise the solver wouldn't work.

 July 4, 2016, 08:00 #19 Senior Member   Join Date: Jan 2013 Posts: 368 Rep Power: 10 Dear All, In twoPhaseEulerFoam solver or many other solvers in OF, we have the terms defined by fvOptions. If in the case system folder, we do not put the fvOptions file, does it mean that all the fvOptions related terms in the source files will not be used or deactived automatically? I am not sure about this and a little how the fvOptions terms work when OpenFOAM is running. Thank you so much.

 July 4, 2016, 08:24 #20 Senior Member   Gerhard Holzinger Join Date: Feb 2012 Location: Austria Posts: 287 Rep Power: 23 You find an answer to your question by playing around. There is a tutorial case of twoPhaseEulerFoam, named injection, which uses the fvOption framework. If you run the tutorial as is, then the solver tells you about the active fvOptions models. See the relevant lines of solver output below. Code: ```No MRF models present Creating finite volume options from "constant/fvOptions" Selecting finite volume options model type scalarSemiImplicitSource Source: massSource1 - selecting cells using points - selected 1 cell(s) with volume 8e-06 Selecting finite volume options model type vectorSemiImplicitSource Source: momentumSource1 - selecting cells using points - selected 1 cell(s) with volume 8e-06 Selecting finite volume options model type scalarSemiImplicitSource Source: energySource1 - selecting cells using points - selected 1 cell(s) with volume 8e-06 Courant Number mean: 0 max: 0 Max Ur Courant Number = 0 Calculating field DDtU1 and DDtU2 Starting time loop``` Then, I renamed the file fvOptions, any name other than fvOptions does the trick. With a non-present fvOptions file, the solver output looked like this: Code: ```No MRF models present No finite volume options present Courant Number mean: 0 max: 0 Max Ur Courant Number = 0 Calculating field DDtU1 and DDtU2 Starting time loop``` openfoammaofnepo likes this.