An eulerian-eulerian approach for oil&gas separator design conference Paper

Turbulence closure 
In the present work, the multiphase turbulence is handled under the two-equation mixture 
𝑘 − 𝜀
model. The governing equations are: 
𝜕
𝜕𝑡
(𝜌
𝑚
𝑘
𝑚
) + ∇ ⋅ (𝜌
𝑚
𝒗
̅
𝑚
𝑘
𝑚
) = ∇ ⋅ [
𝜇
𝑚
𝑡
𝜎
𝑚
∇𝑘
𝑚
] + 𝐺
𝑚
− 𝜌
𝑚
𝜀
𝑚
+ 𝑆
𝑘
𝑚
𝜕
𝜕𝑡
(𝜌
𝑚
𝜀
𝑚
) + ∇ ⋅ (𝜌
𝑚
𝒗
̅
𝑚
𝜀
𝑚
) = ∇ ⋅ [
𝜇
𝑚
𝑡
𝜎
𝑚
∇𝜀
𝑚
] +
𝜀
𝑚
𝑘
𝑚
(𝐶
1
𝐺
𝑚
− 𝐶
2
𝜌
𝑚
𝜀
𝑚
) + 𝐶
3
𝜀
𝑚
𝑘
𝑚
𝑆
𝑘
𝑚
The derivation of this model is based on the summation of a set of 
𝑘 − 𝜀
models written for the single 
phase. In order to relate the turbulent quantities of all phases with respect to one of them, the concept 
of the response coefficient is introduced. Rigorously, the response coefficient is defined as the ratio 
the dispersed phase velocity fluctuations to those of the continuous phase, and is defined as: 
Figure 2: CdRe as function of particle Reynolds number for three 
formulation (Original, Jump and Polynomial)
 

---

6 
𝐶
𝑡,𝑘𝑗
=
𝒗
̅
𝑗
′
𝒗
̅
𝑘
′
As the phase fraction increases, the dominance of the continuous phase on turbulence gradually 
diminishes as that phase becomes confined to thin interstices between dispersed phase elements. 
In the limit of unity dispersed phase fraction, this phase becomes continuous and its turbulence 
becomes the dominant/sole factor. Moreover, experiments [2-3-4] have proved at as the phase 
fraction increases beyond a certain limit, which could be as small as 6 %, the ratio of the dispersed 
to continuous phase fluctuations, 
𝐶
𝑡,𝑘𝑗
,
approaches a constant value close to unity. Based on the 
scaling approach, the turbulent quantities 
𝑘
2
and 
𝑘
3
can be computed from 
𝑘
1
. The same approach 
can be applied to the dissipation energy equation. 
The last aspect to deal with consists in finding a suitable closure for the turbulent kinetic energy 
source term 
𝑆
𝑘
𝑚
. Issa, Behazadi and Rusche [2] derived a model for two-phase flow at high-phase 
fraction: 
𝑆
𝑘,𝑘𝑗
= −𝐴
𝑑,𝑘𝑗
[2𝛼
𝑗
(𝐶
𝑡,𝑘𝑗
− 1)
2
𝑘
𝑘
+ 𝜂
𝑘
∇𝛼
𝑗
⋅ 𝒗
̅
𝑟
]
The previous equation is derived under the assumption that the only contribution to the generation 
of turbulent kinetic energy is represented from the drag force exchange across the interface. Based 
on this idea, our proposal is based on applying the superimposition principle for the drag turbulent 
contribution, exploiting the idea of binary momentum interaction between phases. In this way, the 
source term can be re-written as: 
𝑆
𝑘
= − ∑ 𝐴
𝑑,𝑘𝑗
[2𝛼
𝑗
(𝐶
𝑡,𝑘𝑗
− 1)
2
𝑘
𝑘
+ 𝜂
𝑘
t
∇𝛼
𝑗
⋅ (𝒗
̅
𝑗
− 𝒗
̅
𝑘
)]
3
𝑘=1,𝑗≠𝑘
≈ − ∑ 𝐴
𝑑,𝑘𝑗
[2𝛼
𝑗
(𝐶
𝑡,𝑘𝑗
− 1)
2
𝑘
𝑘
]
3
𝑘=1,𝑗≠𝑘
One of the numerical problem that may arises from this formulation is concerned with the very high 
value that term 
∇𝛼
𝑗
may assume across the interface, leading to serious numerical instabilities. For 
this reason, the gradient of the phase fraction is neglected in the present simulations. This model is 
not claimed to be completely rigorous and is only intended to serve as a preliminary vehicle to test 
the numerical solution procedure until a better closure model becomes available.