Interfacial Composition Model formulation
Summary
The problem is related to the phase equilibrium formulation in reactingEulerFoam. The phase equilibrium for saturated and NRTL models are defined based on the mole-fraction while the solver uses mass fraction basis for its calculations.
Steps to reproduce
Example case
in file src/phaseSystemModels/reactingEulerFoam/interfacialCompositionModels/interfaceCompositionModels/NonRandomTwoLiquid/NonRandomTwoLiquid.C
Line 208 or 215, the code uses the following equation to calculate specie mass fraction in the vapor phase (Yf) : return this->otherThermo_.composition().Y(speciesName) *speciesModel1_->Yf(speciesName, Tf) *gamma1_;
the correct equilibrium based on activity model is based on the mole fraction in pahses:
gamma_i * x_i Ps = y_i * P (Eq. 1)
where both y_i and x_i are mole fraction of component i in liquid and vapor phase. In the code, the mole fraction of the vapor phase is converted to the mass fraction (the return value of Yf method) when we use speciesModel1_->Yf(speciesName, Tf), since in that method the return value is multiplied by Wi/W. However, the conversion of liquid phase mass fraction (otherThermo_.composition().Y(speciesName)) to mole fraction dose not occur. the correct form of the return value is
return this->otherThermo_.composition().Y(speciesName) * otherThermo_.composition().W / otherThermo_.composition().Wi(species1Index_) *speciesModel1_->Yf(speciesName, Tf) *gamma1_;
What is the current bug behaviour?
What is the expected correct behavior?
Relevant logs and/or images
Environment information
Providing details of your set-up can help us identify any issues, e.g. OpenFOAM version : v2006 | V2012 etc Operating system : ubuntu Hardware info : AMD threadripper Compiler : gcc
Possible fixes
in file src/phaseSystemModels/reactingEulerFoam/interfacialCompositionModels/interfaceCompositionModels/NonRandomTwoLiquid/NonRandomTwoLiquid.C
Line 208 should be changed to
dimensionedScalar Wi ( "Wi", dimMass/dimMoles, otherThermo_.composition().Wi(species1Index_) );
return this->otherThermo_.composition().Y(speciesName) * otherThermo_.composition().W /Wi *speciesModel1_->Yf(speciesName, Tf) *gamma1_;
Line 215 should be changed to dimensionedScalar Wi ( "Wi", dimMass/dimMoles, otherThermo_.composition().Wi(species2Index_) ); return this->otherThermo_.composition().Y(speciesName) * otherThermo_.W() / Wi *speciesModel2_->Yf(speciesName, Tf) *gamma2_;